diff --git a/Common/include/code_config.hpp b/Common/include/code_config.hpp index 37b752206b12..3b758a40cf54 100644 --- a/Common/include/code_config.hpp +++ b/Common/include/code_config.hpp @@ -67,6 +67,17 @@ template struct su2conditional { using type = F; template using su2conditional_t = typename su2conditional::type; +/*! \brief Static cast "In" to "Out", in debug builds a dynamic cast is used. */ +template +FORCEINLINE Out su2staticcast_p(In ptr) { + static_assert(std::is_pointer::value, "This expects a pointer"); +#ifndef NDEBUG + return static_cast(ptr); +#else + return dynamic_cast(ptr); +#endif +} + /*--- Detect compilation with OpenMP. ---*/ #if defined(_OPENMP) #define HAVE_OMP diff --git a/SU2_CFD/include/limiters/CLimiterDetails.hpp b/SU2_CFD/include/limiters/CLimiterDetails.hpp index d80a3f407146..b17560763bc4 100644 --- a/SU2_CFD/include/limiters/CLimiterDetails.hpp +++ b/SU2_CFD/include/limiters/CLimiterDetails.hpp @@ -64,22 +64,28 @@ struct CLimiterDetails /*! * \brief Common small functions used by limiters. */ -namespace LimiterHelpers +template +struct LimiterHelpers { - inline passivedouble epsilon() {return std::numeric_limits::epsilon();} + FORCEINLINE static Type epsilon() {return std::numeric_limits::epsilon();} - inline su2double venkatFunction(su2double proj, su2double delta, su2double eps2) + FORCEINLINE static Type venkatFunction(const Type& proj, const Type& delta, const Type& eps2) { - su2double y = delta*(delta+proj) + eps2; + Type y = delta*(delta+proj) + eps2; return (y + delta*proj) / (y + 2*proj*proj); } - inline su2double raisedSine(su2double dist) + FORCEINLINE static Type vanAlbadaFunction(const Type& proj, const Type& delta, const Type& eps) { - su2double factor = 0.5*(1.0+dist+sin(PI_NUMBER*dist)/PI_NUMBER); + return delta * (2*proj + delta) / (4*pow(proj, 2) + pow(delta, 2) + eps); + } + + FORCEINLINE static Type raisedSine(const Type& dist) + { + Type factor = 0.5*(1.0+dist+sin(PI_NUMBER*dist)/PI_NUMBER); return max(0.0, min(factor, 1.0)); } -} +}; /*! @@ -94,7 +100,7 @@ struct CLimiterDetails * \brief Set a small epsilon to avoid divisions by 0. */ template - inline void preprocess(Ts&...) {eps2 = LimiterHelpers::epsilon();} + inline void preprocess(Ts&...) {eps2 = LimiterHelpers<>::epsilon();} /*! * \brief No geometric modification for this kind of limiter. @@ -107,7 +113,7 @@ struct CLimiterDetails */ inline su2double limiterFunction(size_t, su2double proj, su2double delta) const { - return LimiterHelpers::venkatFunction(proj, delta, eps2); + return LimiterHelpers<>::venkatFunction(proj, delta, eps2); } }; @@ -130,7 +136,7 @@ struct CLimiterDetails su2double L = config.GetRefElemLength(); su2double K = config.GetVenkat_LimiterCoeff(); su2double eps1 = fabs(L*K); - eps2 = max(eps1*eps1*eps1, LimiterHelpers::epsilon()); + eps2 = max(eps1*eps1*eps1, LimiterHelpers<>::epsilon()); } /*! @@ -144,7 +150,7 @@ struct CLimiterDetails */ inline su2double limiterFunction(size_t, su2double proj, su2double delta) const { - return LimiterHelpers::venkatFunction(proj, delta, eps2); + return LimiterHelpers<>::venkatFunction(proj, delta, eps2); } }; @@ -229,7 +235,7 @@ struct CLimiterDetails for(size_t iVar = varBegin; iVar < varEnd; ++iVar) { su2double range = sharedMax(iVar) - sharedMin(iVar); - eps2(iVar) = max(pow(K*range, 2), LimiterHelpers::epsilon()); + eps2(iVar) = max(pow(K*range, 2), LimiterHelpers<>::epsilon()); } } @@ -245,7 +251,7 @@ struct CLimiterDetails inline su2double limiterFunction(size_t iVar, su2double proj, su2double delta) const { AD::SetPreaccIn(eps2(iVar)); - return LimiterHelpers::venkatFunction(proj, delta, eps2(iVar)); + return LimiterHelpers<>::venkatFunction(proj, delta, eps2(iVar)); } }; @@ -268,7 +274,7 @@ struct CLimiterDetails su2double L = config.GetRefElemLength(); su2double K = config.GetVenkat_LimiterCoeff(); eps1 = fabs(L*K); - eps2 = max(eps1*eps1*eps1, LimiterHelpers::epsilon()); + eps2 = max(eps1*eps1*eps1, LimiterHelpers<>::epsilon()); } /*! @@ -278,7 +284,7 @@ struct CLimiterDetails { AD::SetPreaccIn(geometry.nodes->GetSharpEdge_Distance(iPoint)); su2double dist = geometry.nodes->GetSharpEdge_Distance(iPoint)/(sharpCoeff*eps1)-1.0; - return LimiterHelpers::raisedSine(dist); + return LimiterHelpers<>::raisedSine(dist); } /*! @@ -286,7 +292,7 @@ struct CLimiterDetails */ inline su2double limiterFunction(size_t, su2double proj, su2double delta) const { - return LimiterHelpers::venkatFunction(proj, delta, eps2); + return LimiterHelpers<>::venkatFunction(proj, delta, eps2); } }; @@ -309,7 +315,7 @@ struct CLimiterDetails su2double L = config.GetRefElemLength(); su2double K = config.GetVenkat_LimiterCoeff(); eps1 = fabs(L*K); - eps2 = max(eps1*eps1*eps1, LimiterHelpers::epsilon()); + eps2 = max(eps1*eps1*eps1, LimiterHelpers<>::epsilon()); } /*! @@ -319,7 +325,7 @@ struct CLimiterDetails { AD::SetPreaccIn(geometry.nodes->GetWall_Distance(iPoint)); su2double dist = geometry.nodes->GetWall_Distance(iPoint)/(sharpCoeff*eps1)-1.0; - return LimiterHelpers::raisedSine(dist); + return LimiterHelpers<>::raisedSine(dist); } /*! @@ -327,6 +333,6 @@ struct CLimiterDetails */ inline su2double limiterFunction(size_t, su2double proj, su2double delta) const { - return LimiterHelpers::venkatFunction(proj, delta, eps2); + return LimiterHelpers<>::venkatFunction(proj, delta, eps2); } }; diff --git a/SU2_CFD/include/solvers/CEulerSolver.hpp b/SU2_CFD/include/solvers/CEulerSolver.hpp index 6c1d0c55aba7..d2d195c2a974 100644 --- a/SU2_CFD/include/solvers/CEulerSolver.hpp +++ b/SU2_CFD/include/solvers/CEulerSolver.hpp @@ -684,12 +684,6 @@ class CEulerSolver : public CFVMFlowSolverBaseSetSolution_New(); } - /*! * \brief Update the solution using a Runge-Kutta scheme. * \param[in] geometry - Geometrical definition of the problem. diff --git a/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp b/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp index d365ffe6b27f..28e128245d3a 100644 --- a/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp +++ b/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp @@ -1102,6 +1102,11 @@ class CFVMFlowSolverBase : public CSolver { public: + /*! + * \brief Set the new solution variables to the current solution value for classical RK. + */ + inline void Set_NewSolution() final { nodes->SetSolution_New(); } + /*! * \brief Load a solution from a restart file. * \param[in] geometry - Geometrical definition of the problem. diff --git a/SU2_CFD/include/solvers/CScalarSolver.inl b/SU2_CFD/include/solvers/CScalarSolver.inl index 31b1280bc5c0..c685d71ff709 100644 --- a/SU2_CFD/include/solvers/CScalarSolver.inl +++ b/SU2_CFD/include/solvers/CScalarSolver.inl @@ -27,6 +27,7 @@ #include "../../../Common/include/parallelization/omp_structure.hpp" #include "../../../Common/include/toolboxes/geometry_toolbox.hpp" #include "../../include/solvers/CScalarSolver.hpp" +#include "../../include/variables/CFlowVariable.hpp" template CScalarSolver::CScalarSolver(bool conservative) : CSolver(), Conservative(conservative) {} @@ -92,10 +93,10 @@ void CScalarSolver::Upwind_Residual(CGeometry* geometry, CSolver** const bool limiterFlow = (config->GetKind_SlopeLimit_Flow() != NO_LIMITER) && (config->GetKind_SlopeLimit_Flow() != VAN_ALBADA_EDGE); - CVariable* flowNodes = solver_container[FLOW_SOL]->GetNodes(); + auto* flowNodes = su2staticcast_p(solver_container[FLOW_SOL]->GetNodes()); /*--- Pick one numerics object per thread. ---*/ - CNumerics* numerics = numerics_container[CONV_TERM + omp_get_thread_num() * MAX_TERMS]; + auto* numerics = numerics_container[CONV_TERM + omp_get_thread_num() * MAX_TERMS]; /*--- Static arrays of MUSCL-reconstructed flow primitives and turbulence variables (thread safety). ---*/ su2double solution_i[MAXNVAR] = {0.0}, flowPrimVar_i[MAXNVARFLOW] = {0.0}; diff --git a/SU2_CFD/include/variables/CEulerVariable.hpp b/SU2_CFD/include/variables/CEulerVariable.hpp index 32f722623f16..4aab05de7910 100644 --- a/SU2_CFD/include/variables/CEulerVariable.hpp +++ b/SU2_CFD/include/variables/CEulerVariable.hpp @@ -27,41 +27,29 @@ #pragma once -#include "CVariable.hpp" +#include "CFlowVariable.hpp" /*! * \class CEulerVariable * \brief Class for defining the variables of the compressible Euler solver. + * \note Primitive variables (T, vx, vy, vz, P, rho, h, c) + * \note Gradients and limiters (T, vx, vy, vz, P, rho) * \ingroup Euler_Equations * \author F. Palacios, T. Economon */ -class CEulerVariable : public CVariable { +class CEulerVariable : public CFlowVariable { public: static constexpr size_t MAXNVAR = 12; -protected: - VectorType Velocity2; /*!< \brief Square of the velocity vector. */ - MatrixType HB_Source; /*!< \brief harmonic balance source term. */ + protected: + /*!< \brief Secondary variables (dPdrho_e, dPde_rho, dTdrho_e, dTde_rho, dmudrho_T, dmudT_rho, dktdrho_T, dktdT_rho) + * in compressible (Euler: 2, NS: 8) flows. */ + MatrixType Secondary; + MatrixType WindGust; /*! < \brief Wind gust value */ MatrixType WindGustDer; /*! < \brief Wind gust derivatives value */ - /*--- Primitive variable definition ---*/ - MatrixType Primitive; /*!< \brief Primitive variables (T, vx, vy, vz, P, rho, h, c) in compressible flows. */ - CVectorOfMatrix Gradient_Primitive; /*!< \brief Gradient of the primitive variables (T, vx, vy, vz, P, rho). */ - CVectorOfMatrix& Gradient_Reconstruction; /*!< \brief Reference to the gradient of the primitive variables for MUSCL reconstruction for the convective term */ - CVectorOfMatrix Gradient_Aux; /*!< \brief Auxiliary structure to store a second gradient for reconstruction, if required. */ - MatrixType Limiter_Primitive; /*!< \brief Limiter of the primitive variables (T, vx, vy, vz, P, rho). */ - - /*--- Secondary variable definition ---*/ - MatrixType Secondary; /*!< \brief Secondary variables (dPdrho_e, dPde_rho, dTdrho_e, dTde_rho, dmudrho_T, dmudT_rho, dktdrho_T, dktdT_rho) in compressible (Euler: 2, NS: 8) flows. */ - - MatrixType Solution_New; /*!< \brief New solution container for Classical RK4. */ - - /*--- NS Variables declared here to make it easier to re-use code between compressible and incompressible solvers. ---*/ - MatrixType Vorticity; /*!< \brief Vorticity of the fluid. */ - VectorType StrainMag; /*!< \brief Magnitude of rate of strain tensor. */ - -public: + public: /*! * \brief Constructor of the class. * \param[in] density - Value of the flow density (initialization value). @@ -73,92 +61,7 @@ class CEulerVariable : public CVariable { * \param[in] config - Definition of the particular problem. */ CEulerVariable(su2double density, const su2double *velocity, su2double energy, - unsigned long npoint, unsigned long ndim, unsigned long nvar, CConfig *config); - - /*! - * \brief Destructor of the class. - */ - ~CEulerVariable() override = default; - - /*! - * \brief Get the new solution of the problem (Classical RK4). - * \param[in] iVar - Index of the variable. - * \return Pointer to the old solution vector. - */ - inline su2double GetSolution_New(unsigned long iPoint, unsigned long iVar) const final { return Solution_New(iPoint,iVar); } - - /*! - * \brief Set the new solution container for Classical RK4. - */ - void SetSolution_New() final; - - /*! - * \brief Add a value to the new solution container for Classical RK4. - * \param[in] iVar - Number of the variable. - * \param[in] val_solution - Value that we want to add to the solution. - */ - inline void AddSolution_New(unsigned long iPoint, unsigned long iVar, su2double val_solution) final { - Solution_New(iPoint,iVar) += val_solution; - } - - /*! - * \brief Get the value of the primitive variables gradient. - * \param[in] iVar - Index of the variable. - * \param[in] iDim - Index of the dimension. - * \return Value of the primitive variables gradient. - */ - inline su2double GetGradient_Primitive(unsigned long iPoint, unsigned long iVar, unsigned long iDim) const final { - return Gradient_Primitive(iPoint,iVar,iDim); - } - - /*! - * \brief Get the primitive variables limiter. - * \return Primitive variables limiter for the entire domain. - */ - inline MatrixType& GetLimiter_Primitive() final {return Limiter_Primitive; } - inline const MatrixType& GetLimiter_Primitive() const final {return Limiter_Primitive; } - - /*! - * \brief Get the value of the primitive variables gradient. - * \param[in] iVar - Index of the variable. - * \return Value of the primitive variables gradient. - */ - inline su2double GetLimiter_Primitive(unsigned long iPoint, unsigned long iVar) const final {return Limiter_Primitive(iPoint,iVar); } - - /*! - * \brief Get the primitive variable gradients for all points. - * \return Reference to primitive variable gradient. - */ - inline CVectorOfMatrix& GetGradient_Primitive() final { return Gradient_Primitive; } - inline const CVectorOfMatrix& GetGradient_Primitive() const final { return Gradient_Primitive; } - - /*! - * \brief Get the reconstruction gradient for primitive variable at all points. - * \return Reference to variable reconstruction gradient. - */ - inline CVectorOfMatrix& GetGradient_Reconstruction() final { return Gradient_Reconstruction; } - inline const CVectorOfMatrix& GetGradient_Reconstruction() const final { return Gradient_Reconstruction; } - - /*! - * \brief Get the value of the primitive variables gradient. - * \return Value of the primitive variables gradient. - */ - inline CMatrixView GetGradient_Primitive(unsigned long iPoint, unsigned long iVar=0) final { - return Gradient_Primitive(iPoint,iVar); - } - - /*! - * \brief Get the value of the primitive variables gradient. - * \return Value of the primitive variables gradient. - */ - inline su2double *GetLimiter_Primitive(unsigned long iPoint) final { return Limiter_Primitive[iPoint]; } - - /*! - * \brief Get the array of the reconstruction variables gradient at a node. - * \param[in] iPoint - Index of the current node. - * \return Array of the reconstruction variables gradient at a node. - */ - inline CMatrixView GetGradient_Reconstruction(unsigned long iPoint) final { return Gradient_Reconstruction[iPoint]; } + unsigned long npoint, unsigned long ndim, unsigned long nvar, const CConfig *config); /*! * \brief A virtual member. @@ -207,43 +110,6 @@ class CEulerVariable : public CVariable { */ void SetSecondaryVar(unsigned long iPoint, CFluidModel *FluidModel) override; - /*! - * \brief Get the primitive variables for all points. - * \return Reference to primitives. - */ - inline const MatrixType& GetPrimitive() const final { return Primitive; } - - /*! - * \brief Get the primitive variables. - * \param[in] iVar - Index of the variable. - * \return Value of the primitive variable for the index iVar. - */ - inline su2double GetPrimitive(unsigned long iPoint, unsigned long iVar) const final { return Primitive(iPoint,iVar); } - - /*! - * \brief Set the value of the primitive variables. - * \param[in] iVar - Index of the variable. - * \param[in] iVar - Index of the variable. - * \return Set the value of the primitive variable for the index iVar. - */ - inline void SetPrimitive(unsigned long iPoint, unsigned long iVar, su2double val_prim) final { Primitive(iPoint,iVar) = val_prim; } - - /*! - * \brief Set the value of the primitive variables. - * \param[in] val_prim - Primitive variables. - * \return Set the value of the primitive variable for the index iVar. - */ - inline void SetPrimitive(unsigned long iPoint, const su2double *val_prim) final { - for (unsigned long iVar = 0; iVar < nPrimVar; iVar++) - Primitive(iPoint,iVar) = val_prim[iVar]; - } - - /*! - * \brief Get the primitive variables of the problem. - * \return Pointer to the primitive variable vector. - */ - inline su2double *GetPrimitive(unsigned long iPoint) final {return Primitive[iPoint]; } - /*! * \brief Get all the secondary variables. */ @@ -254,7 +120,7 @@ class CEulerVariable : public CVariable { * \param[in] iVar - Index of the variable. * \return Value of the secondary variable for the index iVar. */ - inline su2double GetSecondary(unsigned long iPoint, unsigned long iVar) const final {return Secondary(iPoint,iVar); } + inline su2double GetSecondary(unsigned long iPoint, unsigned long iVar) const final { return Secondary(iPoint,iVar); } /*! * \brief Set the value of the secondary variables. @@ -262,7 +128,9 @@ class CEulerVariable : public CVariable { * \param[in] iVar - Index of the variable. * \return Set the value of the secondary variable for the index iVar. */ - inline void SetSecondary(unsigned long iPoint, unsigned long iVar, su2double val_secondary) final {Secondary(iPoint,iVar) = val_secondary; } + inline void SetSecondary(unsigned long iPoint, unsigned long iVar, su2double val_secondary) final { + Secondary(iPoint,iVar) = val_secondary; + } /*! * \brief Set the value of the secondary variables. @@ -297,12 +165,6 @@ class CEulerVariable : public CVariable { return temperature <= 0.0; } - /*! - * \brief Get the norm 2 of the velocity. - * \return Norm 2 of the velocity vector. - */ - inline su2double GetVelocity2(unsigned long iPoint) const final { return Velocity2(iPoint); } - /*! * \brief Get the flow pressure. * \return Value of the flow pressure. @@ -387,22 +249,6 @@ class CEulerVariable : public CVariable { for (unsigned long iDim = 0; iDim < nDim; iDim++) Res_TruncError(iPoint,iDim+1) = 0.0; } - /*! - * \brief Set the harmonic balance source term. - * \param[in] iVar - Index of the variable. - * \param[in] val_solution - Value of the harmonic balance source term. for the index iVar. - */ - inline void SetHarmonicBalance_Source(unsigned long iPoint, unsigned long iVar, su2double val_source) final { - HB_Source(iPoint,iVar) = val_source; - } - - /*! - * \brief Get the harmonic balance source term. - * \param[in] iVar - Index of the variable. - * \return Value of the harmonic balance source term for the index iVar. - */ - inline su2double GetHarmonicBalance_Source(unsigned long iPoint, unsigned long iVar) const final { return HB_Source(iPoint,iVar); } - /*! * \brief Get the value of the wind gust * \return Value of the wind gust @@ -433,19 +279,6 @@ class CEulerVariable : public CVariable { WindGustDer(iPoint,iDim) = val_WindGustDer[iDim]; } - /*! - * \brief Get the value of the vorticity. - * \return Value of the vorticity. - */ - inline su2double *GetVorticity(unsigned long iPoint) final { return Vorticity[iPoint]; } - - /*! - * \brief Get the value of the magnitude of rate of strain. - * \return Value of the rate of strain magnitude. - */ - inline su2double GetStrainMag(unsigned long iPoint) const final { return StrainMag(iPoint); } - inline su2activevector& GetStrainMag() { return StrainMag; } - /*! * \brief Specify a vector to set the velocity components of the solution. Multiplied by density for compressible cases. * \param[in] iPoint - Point index. diff --git a/SU2_CFD/include/variables/CFlowVariable.hpp b/SU2_CFD/include/variables/CFlowVariable.hpp new file mode 100644 index 000000000000..f65d36087381 --- /dev/null +++ b/SU2_CFD/include/variables/CFlowVariable.hpp @@ -0,0 +1,252 @@ +/*! + * \file CFlowVariable.hpp + * \brief Class for defining the common variables of flow solvers. + * \version 7.2.0 "Blackbird" + * + * SU2 Project Website: https://su2code.github.io + * + * The SU2 Project is maintained by the SU2 Foundation + * (http://su2foundation.org) + * + * Copyright 2012-2021, 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 . + */ + +#pragma once + +#include "CVariable.hpp" + +/*! + * \class CFlowVariable + * \brief Class for defining the common variables of flow solvers. + */ +class CFlowVariable : public CVariable { + protected: + /*--- Primitive variable definition. ---*/ + MatrixType Primitive; /*!< \brief Primitive variables. */ + CVectorOfMatrix Gradient_Primitive; /*!< \brief Gradient of the primitive variables. */ + CVectorOfMatrix& Gradient_Reconstruction; /*!< \brief Reference to the gradient of the primitive variables for MUSCL + reconstruction for the convective term */ + CVectorOfMatrix + Gradient_Aux; /*!< \brief Auxiliary structure to store a second gradient for reconstruction, if required. */ + MatrixType Limiter_Primitive; /*!< \brief Limiter of the primitive variables. */ + VectorType Velocity2; /*!< \brief Squared norm of velocity. */ + + MatrixType Solution_New; /*!< \brief New solution container for Classical RK4. */ + MatrixType HB_Source; /*!< \brief harmonic balance source term. */ + + /*--- NS Variables declared here to make it easier to re-use code between compressible and incompressible solvers. + * ---*/ + MatrixType Vorticity; /*!< \brief Vorticity of the flow field. */ + VectorType StrainMag; /*!< \brief Magnitude of rate of strain tensor. */ + + /*! + * \brief Constructor of the class. + * \note This class is not meant to be instantiated directly, it is only a building block. + * \param[in] npoint - Number of points/nodes/vertices in the domain. + * \param[in] ndim - Number of dimensions. + * \param[in] nvar - Number of variables. + * \param[in] nprimvar - Number of primitive variables. + * \param[in] nprimvargrad - Number of primitive variables for which to compute gradients. + * \param[in] config - Definition of the particular problem. + */ + CFlowVariable(unsigned long npoint, unsigned long ndim, unsigned long nvar, unsigned long nprimvar, + unsigned long nprimvargrad, const CConfig* config); + + public: + /*! + * \brief Get a primitive variable. + * \param[in] iPoint - Point index. + * \param[in] iVar - Index of the variable. + * \return Value of the primitive variable for the index iVar. + */ + inline su2double GetPrimitive(unsigned long iPoint, unsigned long iVar) const final { + return Primitive(iPoint, iVar); + } + + /*! + * \brief Get all the primitive variables of the problem. + * \param[in] iPoint - Point index. + * \return Pointer to the primitive variable vector. + */ + inline su2double* GetPrimitive(unsigned long iPoint) final { return Primitive[iPoint]; } + + /*! + * \brief Get the primitive variables for all points. + * \return Reference to primitives. + */ + inline const MatrixType& GetPrimitive() const final { return Primitive; } + + /*! + * \brief Set the value of the primitive variables. + * \param[in] iPoint - Point index. + * \param[in] iVar - Index of the variable. + * \param[in] val_prim - Value of the variable. + * \return Set the value of the primitive variable for the index iVar. + */ + inline void SetPrimitive(unsigned long iPoint, unsigned long iVar, su2double val_prim) final { + Primitive(iPoint, iVar) = val_prim; + } + + /*! + * \brief Set the value of the primitive variables. + * \param[in] iPoint - Point index. + * \param[in] val_prim - Values of primitive variables. + * \return Set the value of the primitive variable for the index iVar. + */ + inline void SetPrimitive(unsigned long iPoint, const su2double* val_prim) final { + for (unsigned long iVar = 0; iVar < nPrimVar; iVar++) Primitive(iPoint, iVar) = val_prim[iVar]; + } + + /*! + * \brief Get the squared norm of the velocity. + * \param[in] iPoint - Point index. + * \return Squared norm of the velocity vector. + */ + inline su2double GetVelocity2(unsigned long iPoint) const final { return Velocity2(iPoint); } + + /*! + * \brief Get the value of the primitive variables gradient. + * \param[in] iPoint - Point index. + * \param[in] iVar - Index of the variable. + * \param[in] iDim - Index of the dimension. + * \return Value of the primitive variables gradient. + */ + inline su2double GetGradient_Primitive(unsigned long iPoint, unsigned long iVar, unsigned long iDim) const final { + return Gradient_Primitive(iPoint, iVar, iDim); + } + + /*! + * \brief Get the value of the primitive variables gradient. + * \param[in] iPoint - Point index. + * \param[in] iVarStart - Offset from which to start reading. + * \return Value of the primitive variables gradient. + */ + inline CMatrixView GetGradient_Primitive(unsigned long iPoint, unsigned long iVarStart = 0) final { + return Gradient_Primitive(iPoint, iVarStart); + } + + /*! + * \brief Get the primitive variable gradients for all points. + * \return Reference to primitive variable gradient. + */ + inline CVectorOfMatrix& GetGradient_Primitive() final { return Gradient_Primitive; } + inline const CVectorOfMatrix& GetGradient_Primitive() const final { return Gradient_Primitive; } + + /*! + * \brief Get the array of the reconstruction variables gradient at a node. + * \param[in] iPoint - Point index. + * \return Array of the reconstruction variables gradient at a node. + */ + inline CMatrixView GetGradient_Reconstruction(unsigned long iPoint) final { + return Gradient_Reconstruction[iPoint]; + } + + /*! + * \brief Get the reconstruction gradient for primitive variable at all points. + * \return Reference to variable reconstruction gradient. + */ + inline CVectorOfMatrix& GetGradient_Reconstruction() final { return Gradient_Reconstruction; } + inline const CVectorOfMatrix& GetGradient_Reconstruction() const final { return Gradient_Reconstruction; } + + /*! + * \brief Get the value of the primitive variables gradient. + * \param[in] iPoint - Point index. + * \param[in] iVar - Index of the variable. + * \return Value of the primitive variables gradient. + */ + inline su2double GetLimiter_Primitive(unsigned long iPoint, unsigned long iVar) const final { + return Limiter_Primitive(iPoint, iVar); + } + + /*! + * \brief Get the value of the primitive variables gradient. + * \param[in] iPoint - Point index. + * \return Value of the primitive variables gradient. + */ + inline su2double* GetLimiter_Primitive(unsigned long iPoint) final { return Limiter_Primitive[iPoint]; } + + /*! + * \brief Get the primitive variables limiter. + * \return Primitive variables limiter for the entire domain. + */ + inline MatrixType& GetLimiter_Primitive() final { return Limiter_Primitive; } + inline const MatrixType& GetLimiter_Primitive() const final { return Limiter_Primitive; } + + /*! + * \brief Get the new solution of the problem (Classical RK4). + * \param[in] iPoint - Point index. + * \param[in] iVar - Index of the variable. + * \return Pointer to the old solution vector. + */ + inline su2double GetSolution_New(unsigned long iPoint, unsigned long iVar) const final { + return Solution_New(iPoint, iVar); + } + + /*! + * \brief Add a value to the new solution container for Classical RK4. + * \param[in] iPoint - Point index. + * \param[in] iVar - Index of the variable. + * \param[in] val_solution - Value that we want to add to the new solution. + */ + inline void AddSolution_New(unsigned long iPoint, unsigned long iVar, su2double val_solution) final { + Solution_New(iPoint, iVar) += val_solution; + } + + /*! + * \brief Set the new solution container for Classical RK4. + */ + void SetSolution_New() final; + + /*! + * \brief Get the harmonic balance source term. + * \param[in] iPoint - Point index. + * \param[in] iVar - Index of the variable. + * \return Value of the harmonic balance source term for the index iVar. + */ + inline su2double GetHarmonicBalance_Source(unsigned long iPoint, unsigned long iVar) const final { + return HB_Source(iPoint, iVar); + } + + /*! + * \brief Set the harmonic balance source term. + * \param[in] iPoint - Point index. + * \param[in] iVar - Index of the variable. + * \param[in] val_solution - Value of the harmonic balance source term. for the index iVar. + */ + inline void SetHarmonicBalance_Source(unsigned long iPoint, unsigned long iVar, su2double val_source) final { + HB_Source(iPoint, iVar) = val_source; + } + + /*! + * \brief Get the values of the vorticity (3 values also in 2D). + * \param[in] iPoint - Point index. + * \return Vorticity array. + */ + inline su2double* GetVorticity(unsigned long iPoint) final { return Vorticity[iPoint]; } + + /*! + * \brief Get the magnitude of rate of strain. + * \param[in] iPoint - Point index. + * \return Value of magnitude. + */ + inline su2double GetStrainMag(unsigned long iPoint) const final { return StrainMag(iPoint); } + + /*! + * \brief Get the entire vector of the rate of strain magnitude. + * \return Vector of magnitudes. + */ + inline su2activevector& GetStrainMag() { return StrainMag; } +}; diff --git a/SU2_CFD/include/variables/CIncEulerVariable.hpp b/SU2_CFD/include/variables/CIncEulerVariable.hpp index 68b984650d7b..0d75dbe36868 100644 --- a/SU2_CFD/include/variables/CIncEulerVariable.hpp +++ b/SU2_CFD/include/variables/CIncEulerVariable.hpp @@ -27,34 +27,24 @@ #pragma once -#include "CVariable.hpp" +#include "CFlowVariable.hpp" /*! * \class CIncEulerVariable * \brief Class for defining the variables of the incompressible Euler solver. + * \note Primitive variables (P, vx, vy, vz, T, rho, beta, lamMu, EddyMu, Kt_eff, Cp, Cv) + * \note Gradients of primitives (P, vx, vy, vz, T, rho, beta) * \ingroup Euler_Equations * \author F. Palacios, T. Economon, T. Albring */ -class CIncEulerVariable : public CVariable { +class CIncEulerVariable : public CFlowVariable { public: static constexpr size_t MAXNVAR = 12; -protected: - VectorType Velocity2; /*!< \brief Square of the velocity vector. */ - MatrixType Primitive; /*!< \brief Primitive variables (P, vx, vy, vz, T, rho, beta, lamMu, EddyMu, Kt_eff, Cp, Cv) in incompressible flows. */ - CVectorOfMatrix Gradient_Primitive; /*!< \brief Gradient of the primitive variables (P, vx, vy, vz, T, rho, beta). */ - CVectorOfMatrix& Gradient_Reconstruction; /*!< \brief Reference to the gradient of the primitive variables for MUSCL reconstruction for the convective term */ - CVectorOfMatrix Gradient_Aux; /*!< \brief Auxiliary structure to store a second gradient for reconstruction, if required. */ - MatrixType Limiter_Primitive; /*!< \brief Limiter of the primitive variables (P, vx, vy, vz, T, rho, beta). */ - - /*--- NS Variables declared here to make it easier to re-use code between compressible and incompressible solvers. ---*/ - MatrixType Vorticity; /*!< \brief Vorticity of the fluid. */ - VectorType StrainMag; /*!< \brief Magnitude of rate of strain tensor. */ - + protected: VectorType Streamwise_Periodic_RecoveredPressure, /*!< \brief Recovered/Physical pressure [Pa] for streamwise periodic flow. */ Streamwise_Periodic_RecoveredTemperature; /*!< \brief Recovered/Physical temperature [K] for streamwise periodic flow. */ - -public: + public: /*! * \brief Constructor of the class. * \param[in] val_pressure - value of the pressure. @@ -66,77 +56,7 @@ class CIncEulerVariable : public CVariable { * \param[in] config - Definition of the particular problem. */ CIncEulerVariable(su2double pressure, const su2double *velocity, su2double temperature, - unsigned long npoint, unsigned long ndim, unsigned long nvar, CConfig *config); - - /*! - * \brief Destructor of the class. - */ - ~CIncEulerVariable() override = default; - - /*! - * \brief Get the primitive variable gradients for all points. - * \return Reference to primitive variable gradient. - */ - inline CVectorOfMatrix& GetGradient_Primitive() final { return Gradient_Primitive; } - inline const CVectorOfMatrix& GetGradient_Primitive() const final { return Gradient_Primitive; } - - /*! - * \brief Get the reconstruction gradient for primitive variable at all points. - * \return Reference to variable reconstruction gradient. - */ - inline CVectorOfMatrix& GetGradient_Reconstruction() final { return Gradient_Reconstruction; } - inline const CVectorOfMatrix& GetGradient_Reconstruction() const final { return Gradient_Reconstruction; } - - /*! - * \brief Get the value of the primitive variables gradient. - * \param[in] iPoint - Point index. - * \param[in] iVar - Index of the variable. - * \param[in] iDim - Index of the dimension. - * \return Value of the primitive variables gradient. - */ - inline su2double GetGradient_Primitive(unsigned long iPoint, unsigned long iVar, unsigned long iDim) const final { - return Gradient_Primitive(iPoint,iVar,iDim); - } - - /*! - * \brief Get the primitive variables limiter. - * \return Primitive variables limiter for the entire domain. - */ - inline MatrixType& GetLimiter_Primitive() final {return Limiter_Primitive; } - inline const MatrixType& GetLimiter_Primitive() const final {return Limiter_Primitive; } - - /*! - * \brief Get the value of the primitive variables gradient. - * \param[in] iPoint - Point index. - * \param[in] iVar - Index of the variable. - * \return Value of the primitive variables gradient. - */ - inline su2double GetLimiter_Primitive(unsigned long iPoint, unsigned long iVar) const final { - return Limiter_Primitive(iPoint,iVar); - } - - /*! - * \brief Get the value of the primitive variables gradient. - * \param[in] iPoint - Point index. - * \return Value of the primitive variables gradient. - */ - inline CMatrixView GetGradient_Primitive(unsigned long iPoint, unsigned long iVar=0) final { - return Gradient_Primitive(iPoint,iVar); - } - - /*! - * \brief Get the value of the primitive variables gradient. - * \param[in] iPoint - Point index. - * \return Value of the primitive variables gradient. - */ - inline su2double *GetLimiter_Primitive(unsigned long iPoint) final { return Limiter_Primitive[iPoint]; } - - /*! - * \brief Get the array of the reconstruction variables gradient at a node. - * \param[in] iPoint - Index of the current node. - * \return Array of the reconstruction variables gradient at a node. - */ - inline CMatrixView GetGradient_Reconstruction(unsigned long iPoint) final { return Gradient_Reconstruction[iPoint]; } + unsigned long npoint, unsigned long ndim, unsigned long nvar, const CConfig *config); /*! * \brief Set the value of the pressure. @@ -144,46 +64,6 @@ class CIncEulerVariable : public CVariable { */ inline void SetPressure(unsigned long iPoint) final { Primitive(iPoint,0) = Solution(iPoint,0); } - /*! - * \brief Get the primitive variables for all points. - * \return Reference to primitives. - */ - inline const MatrixType& GetPrimitive() const final { return Primitive; } - - /*! - * \brief Get the primitive variables. - * \param[in] iPoint - Point index. - * \param[in] iVar - Index of the variable. - * \return Value of the primitive variable for the index iVar. - */ - inline su2double GetPrimitive(unsigned long iPoint, unsigned long iVar) const final { return Primitive(iPoint,iVar); } - - /*! - * \brief Set the value of the primitive variables. - * \param[in] iPoint - Point index. - * \param[in] iVar - Index of the variable. - * \param[in] iVar - Index of the variable. - * \return Set the value of the primitive variable for the index iVar. - */ - inline void SetPrimitive(unsigned long iPoint, unsigned long iVar, su2double val_prim) final { Primitive(iPoint,iVar) = val_prim; } - - /*! - * \brief Set the value of the primitive variables. - * \param[in] iPoint - Point index. - * \param[in] val_prim - Primitive variables. - * \return Set the value of the primitive variable for the index iVar. - */ - inline void SetPrimitive(unsigned long iPoint, const su2double *val_prim) final { - for (unsigned long iVar = 0; iVar < nPrimVar; iVar++) Primitive(iPoint,iVar) = val_prim[iVar]; - } - - /*! - * \brief Get the primitive variables of the problem. - * \param[in] iPoint - Point index. - * \return Pointer to the primitive variable vector. - */ - inline su2double *GetPrimitive(unsigned long iPoint) final { return Primitive[iPoint]; } - /*! * \brief Set the value of the density for the incompressible flows. * \param[in] iPoint - Point index. @@ -222,12 +102,6 @@ class CIncEulerVariable : public CVariable { */ inline void SetBetaInc2(unsigned long iPoint, su2double val_betainc2) final { Primitive(iPoint,nDim+3) = val_betainc2; } - /*! - * \brief Get the norm 2 of the velocity. - * \return Norm 2 of the velocity vector. - */ - inline su2double GetVelocity2(unsigned long iPoint) const final { return Velocity2(iPoint); } - /*! * \brief Get the flow pressure. * \return Value of the flow pressure. @@ -315,19 +189,6 @@ class CIncEulerVariable : public CVariable { */ inline su2double GetSpecificHeatCv(unsigned long iPoint) const final { return Primitive(iPoint, nDim+8); } - /*! - * \brief Get the value of the vorticity. - * \return Value of the vorticity. - */ - inline su2double *GetVorticity(unsigned long iPoint) final { return Vorticity[iPoint]; } - - /*! - * \brief Get the value of the magnitude of rate of strain. - * \return Value of the rate of strain magnitude. - */ - inline su2double GetStrainMag(unsigned long iPoint) const final { return StrainMag(iPoint); } - inline su2activevector& GetStrainMag() { return StrainMag; } - /*! * \brief Set the recovered pressure for streamwise periodic flow. * \param[in] iPoint - Point index. diff --git a/SU2_CFD/include/variables/CIncNSVariable.hpp b/SU2_CFD/include/variables/CIncNSVariable.hpp index 0ffedfeb5a9c..4245f64c44f6 100644 --- a/SU2_CFD/include/variables/CIncNSVariable.hpp +++ b/SU2_CFD/include/variables/CIncNSVariable.hpp @@ -52,12 +52,7 @@ class CIncNSVariable final : public CIncEulerVariable { * \param[in] config - Definition of the particular problem. */ CIncNSVariable(su2double pressure, const su2double *velocity, su2double temperature, - unsigned long npoint, unsigned long ndim, unsigned long nvar, CConfig *config); - - /*! - * \brief Destructor of the class. - */ - ~CIncNSVariable() override = default; + unsigned long npoint, unsigned long ndim, unsigned long nvar, const CConfig *config); /*! * \brief Set the laminar viscosity. diff --git a/SU2_CFD/include/variables/CNEMOEulerVariable.hpp b/SU2_CFD/include/variables/CNEMOEulerVariable.hpp index fd728c3e7eeb..9f5f65a9c17d 100644 --- a/SU2_CFD/include/variables/CNEMOEulerVariable.hpp +++ b/SU2_CFD/include/variables/CNEMOEulerVariable.hpp @@ -27,48 +27,34 @@ #pragma once -#include "CVariable.hpp" +#include "CFlowVariable.hpp" #include "../fluid/CNEMOGas.hpp" #include "../../../Common/include/toolboxes/geometry_toolbox.hpp" /*! * \class CNEMOEulerVariable * \brief Main class for defining the variables of the NEMO Euler's solver. + * \note Primitive variables (rhos_s, T, Tve, vx, vy, vw, P, rho, h, rhoCvtr, rhoCvve) * \ingroup Euler_Equations * \author S. R. Copeland, F. Palacios, W. Maier, C. Garbacz */ -class CNEMOEulerVariable : public CVariable { -public: +class CNEMOEulerVariable : public CFlowVariable { + public: static constexpr size_t MAXNVAR = 25; -protected: - + protected: bool ionization; /*!< \brief Presence of charged species in gas mixture. */ bool monoatomic = false; /*!< \brief Presence of single species gas. */ - VectorType Velocity2; /*!< \brief Square of the velocity vector. */ MatrixType Precond_Beta; /*!< \brief Low Mach number preconditioner value, Beta. */ - CVectorOfMatrix& Gradient_Reconstruction; /*!< \brief Reference to the gradient of the conservative variables for MUSCL reconstruction for the convective term */ - CVectorOfMatrix Gradient_Aux; /*!< \brief Auxiliary structure to store a second gradient for reconstruction, if required. */ - /*--- Primitive variable definition ---*/ - MatrixType Primitive; /*!< \brief Primitive variables (rhos_s, T, Tve, ...) in compressible flows. */ MatrixType Primitive_Aux; /*!< \brief Primitive auxiliary variables (Y_s, T, Tve, ...) in compressible flows. */ - CVectorOfMatrix Gradient_Primitive; /*!< \brief Gradient of the primitive variables (rhos_s, T, Tve, ...). */ - MatrixType Limiter_Primitive; /*!< \brief Limiter of the primitive variables (rhos_s, T, Tve, ...). */ /*--- Secondary variable definition ---*/ MatrixType Secondary; /*!< \brief Primitive variables (T, vx, vy, vz, P, rho, h, c) in compressible flows. */ CVectorOfMatrix Gradient_Secondary; /*!< \brief Gradient of the primitive variables (T, vx, vy, vz, P, rho). */ - /*--- New solution container for Classical RK4 ---*/ - MatrixType Solution_New; /*!< \brief New solution container for Classical RK4. */ - - /*--- NS Variables declared here to make it easier to re-use code between compressible and incompressible solvers. ---*/ - MatrixType Vorticity; /*!< \brief Vorticity of the fluid. */ - VectorType StrainMag; /*!< \brief Magnitude of rate of strain tensor. */ - /*--- Other Necessary Variable Definition ---*/ MatrixType dPdU; /*!< \brief Partial derivative of pressure w.r.t. conserved variables. */ MatrixType dTdU; /*!< \brief Partial derivative of temperature w.r.t. conserved variables. */ @@ -87,8 +73,7 @@ class CNEMOEulerVariable : public CVariable { su2double Tve_Freestream; /*!< \brief Freestream vib-el temperature. */ const bool implicit; /*!< \brief Implicit flag. */ -public: - + public: /*! * \brief Constructor of the class. * \param[in] val_pressure - Value of the flow pressure (initialization value). @@ -104,81 +89,16 @@ class CNEMOEulerVariable : public CVariable { * \param[in] config - Definition of the particular problem. */ CNEMOEulerVariable(su2double val_pressure, const su2double *val_massfrac, - su2double *val_mach, su2double val_temperature, + const su2double *val_mach, su2double val_temperature, su2double val_temperature_ve, unsigned long npoint, unsigned long ndim, unsigned long nvar, unsigned long nvalprim, - unsigned long nvarprimgrad, CConfig *config, CNEMOGas *fluidmodel); - - /*! - * \brief Destructor of the class. - */ - ~CNEMOEulerVariable() override = default; + unsigned long nvarprimgrad, const CConfig *config, CNEMOGas *fluidmodel); /*---------------------------------------*/ /*--- U,V,S Routines ---*/ /*---------------------------------------*/ - /*! - * \brief Get the new solution of the problem (Classical RK4). - * \param[in] iVar - Index of the variable. - * \return Pointer to the old solution vector. - */ - inline su2double GetSolution_New(unsigned long iPoint, unsigned long iVar) const final { return Solution_New(iPoint,iVar); } - - /*! - * \brief Set the new solution container for Classical RK4. - */ - void SetSolution_New() final; - - /*! - * \brief Add a value to the new solution container for Classical RK4. - * \param[in] iVar - Number of the variable. - * \param[in] val_solution - Value that we want to add to the solution. - */ - inline void AddSolution_New(unsigned long iPoint, unsigned long iVar, su2double val_solution) final { - Solution_New(iPoint,iVar) += val_solution; - } - - /*! - * \brief Set the value of the primitive variables. - * \param[in] iVar - Index of the variable. - * \param[in] iVar - Index of the variable. - * \return Set the value of the primitive variable for the index iVar. - */ - inline void SetPrimitive(unsigned long iPoint, unsigned long iVar, su2double val_prim) final { Primitive(iPoint,iVar) = val_prim; } - - /*! - * \brief Set the value of the primitive variables. - * \param[in] val_prim - Primitive variables. - * \return Set the value of the primitive variable for the index iVar. - */ - inline void SetPrimitive(unsigned long iPoint, const su2double *val_prim) final { - for (unsigned long iVar = 0; iVar < nPrimVar; iVar++) - Primitive(iPoint,iVar) = val_prim[iVar]; - } - - /*! - * \brief Get the primitive variables limiter. - * \return Primitive variables limiter for the entire domain. - */ - inline MatrixType& GetLimiter_Primitive() final {return Limiter_Primitive; } - inline const MatrixType& GetLimiter_Primitive() const final {return Limiter_Primitive; } - - /*! - * \brief Set the gradient of the primitive variables. - * \param[in] iVar - Index of the variable. - * \param[in] iDim - Index of the dimension. - * \param[in] value - Value of the gradient. - */ - inline su2double GetLimiter_Primitive(unsigned long iPoint, unsigned long iVar) const final {return Limiter_Primitive(iPoint,iVar); } - - /*! - * \brief Get the value of the primitive variables gradient. - * \return Value of the primitive variables gradient. - */ - inline su2double *GetLimiter_Primitive(unsigned long iPoint) final { return Limiter_Primitive[iPoint]; } - /*! * \brief Set the value of the primitive variables. * \param[in] iVar - Index of the variable. @@ -205,32 +125,11 @@ class CNEMOEulerVariable : public CVariable { */ inline void SetPrimitive_Aux(unsigned long iPoint, unsigned long iVar, su2double val_prim) { Primitive_Aux(iPoint,iVar) = val_prim; } - - /*! - * \brief Get the primitive variables. - * \param[in] iVar - Index of the variable. - * \return Value of the primitive variable for the index iVar. - */ - inline su2double GetPrimitive(unsigned long iPoint, unsigned long iVar) const final { return Primitive(iPoint,iVar); } - - /*! - * \brief Get the primitive variables of the problem. - * \return Pointer to the primitive variable vector. - */ - inline su2double *GetPrimitive(unsigned long iPoint) final {return Primitive[iPoint]; } - /*! * \brief Get the primitive variables for all points. * \return Reference to primitives. */ - inline const MatrixType& GetPrimitive() const final { return Primitive; } - - /*! - * \brief Get the primitive variables for all points. - * \return Reference to primitives. - */ - inline const MatrixType& GetPrimitive_Aux(void) const { return Primitive_Aux; } - + inline const MatrixType& GetPrimitive_Aux() const { return Primitive_Aux; } /*! * \brief Get the primitive variables. @@ -245,59 +144,6 @@ class CNEMOEulerVariable : public CVariable { */ inline su2double *GetSecondary(unsigned long iPoint) final { return Secondary[iPoint]; } - /*---------------------------------------*/ - /*--- Gradient Routines ---*/ - /*---------------------------------------*/ - - /*! - * \brief Get the reconstruction gradient for primitive variable at all points. - * \return Reference to variable reconstruction gradient. - */ - inline CVectorOfMatrix& GetGradient_Reconstruction() final { return Gradient_Reconstruction; } - inline const CVectorOfMatrix& GetGradient_Reconstruction() const final { return Gradient_Reconstruction; } - - /*! - * \brief Get the array of the reconstruction variables gradient at a node. - * \param[in] iPoint - Index of the current node. - * \return Array of the reconstruction variables gradient at a node. - */ - inline CMatrixView GetGradient_Reconstruction(unsigned long iPoint) final { return Gradient_Reconstruction[iPoint]; } - - /*! - * \brief Subtract value to the gradient of the primitive variables. - * \param[in] iVar - Index of the variable. - * \param[in] iDim - Index of the dimension. - * \param[in] value - Value to subtract to the gradient of the primitive variables. - */ - inline void SubtractGradient_Primitive(unsigned long iPoint, unsigned long iVar, unsigned long iDim, su2double value) { - Gradient_Primitive(iPoint,iVar,iDim) -= value; - } - - /*! - * \brief Get the value of the primitive variables gradient. - * \param[in] iVar - Index of the variable. - * \param[in] iDim - Index of the dimension. - * \return Value of the primitive variables gradient. - */ - inline su2double GetGradient_Primitive(unsigned long iPoint, unsigned long iVar, unsigned long iDim) const final { - return Gradient_Primitive(iPoint,iVar,iDim); - } - - /*! - * \brief Get the value of the primitive variables gradient. - * \return Value of the primitive variables gradient. - */ - inline CMatrixView GetGradient_Primitive(unsigned long iPoint, unsigned long iVar=0) final { - return Gradient_Primitive(iPoint,iVar); - } - - /*! - * \brief Get the primitive variable gradients for all points. - * \return Reference to primitive variable gradient. - */ - inline CVectorOfMatrix& GetGradient_Primitive() final { return Gradient_Primitive; } - inline const CVectorOfMatrix& GetGradient_Primitive() const final { return Gradient_Primitive; } - /*! * \brief Set all the primitive variables for compressible flows. */ @@ -320,12 +166,6 @@ class CNEMOEulerVariable : public CVariable { */ void SetVelocity2(unsigned long iPoint) final; - /*! - * \brief Get the norm 2 of the velocity. - * \return Norm 2 of the velocity vector. - */ - inline su2double GetVelocity2(unsigned long iPoint) const final { return Velocity2(iPoint); } - /*! * \brief Get the flow pressure. * \return Value of the flow pressure. @@ -483,52 +323,52 @@ class CNEMOEulerVariable : public CVariable { /*! * \brief Retrieves the value of the species density in the primitive variable vector. */ - inline unsigned short GetRhosIndex(void) { return RHOS_INDEX; } + inline unsigned short GetRhosIndex() const { return RHOS_INDEX; } /*! * \brief Retrieves the value of the total density in the primitive variable vector. */ - inline unsigned short GetRhoIndex(void) { return RHO_INDEX; } + inline unsigned short GetRhoIndex() const { return RHO_INDEX; } /*! * \brief Retrieves the value of the pressure in the primitive variable vector. */ - inline unsigned short GetPIndex(void) { return P_INDEX; } + inline unsigned short GetPIndex() const { return P_INDEX; } /*! * \brief Retrieves the value of the in temperature the primitive variable vector. */ - inline unsigned short GetTIndex(void) { return T_INDEX; } + inline unsigned short GetTIndex() const { return T_INDEX; } /*! * \brief Retrieves the value of the vibe-elec temperature in the primitive variable vector. */ - inline unsigned short GetTveIndex(void) { return TVE_INDEX; } + inline unsigned short GetTveIndex() const { return TVE_INDEX; } /*! * \brief Retrieves the value of the velocity in the primitive variable vector. */ - inline unsigned short GetVelIndex(void) { return VEL_INDEX; } + inline unsigned short GetVelIndex() const { return VEL_INDEX; } /*! * \brief Retrieves the value of the enthalpy in the primitive variable vector. */ - inline unsigned short GetHIndex(void) { return H_INDEX; } + inline unsigned short GetHIndex() const { return H_INDEX; } /*! * \brief Retrieves the value of the soundspeed in the primitive variable vector. */ - inline unsigned short GetAIndex(void) { return A_INDEX; } + inline unsigned short GetAIndex() const { return A_INDEX; } /*! * \brief Retrieves the value of the RhoCvtr in the primitive variable vector. */ - inline unsigned short GetRhoCvtrIndex(void) { return RHOCVTR_INDEX; } + inline unsigned short GetRhoCvtrIndex() const { return RHOCVTR_INDEX; } /*! * \brief Retrieves the value of the RhoCvve in the primitive variable vector. */ - inline unsigned short GetRhoCvveIndex(void) { return RHOCVVE_INDEX; } + inline unsigned short GetRhoCvveIndex() const { return RHOCVVE_INDEX; } /*! * \brief Specify a vector to set the velocity components of the solution. Multiplied by density for compressible cases. @@ -540,19 +380,6 @@ class CNEMOEulerVariable : public CVariable { Solution(iPoint, nSpecies+iDim) = Primitive(iPoint,RHO_INDEX) * val_vector[iDim]; } - /*! - * \brief Get the value of the vorticity. - * \return Value of the vorticity. - */ - inline su2double *GetVorticity(unsigned long iPoint) final { return Vorticity[iPoint]; } - - /*! - * \brief Get the value of the magnitude of rate of strain. - * \return Value of the rate of strain magnitude. - */ - inline su2double GetStrainMag(unsigned long iPoint) const final { return StrainMag(iPoint); } - inline su2activevector& GetStrainMag() { return StrainMag; } - /*! * \brief Set the momentum part of the truncation error to zero. * \param[in] iPoint - Point index. diff --git a/SU2_CFD/include/variables/CNEMONSVariable.hpp b/SU2_CFD/include/variables/CNEMONSVariable.hpp index f124700d74e1..5fd0e0f24e12 100644 --- a/SU2_CFD/include/variables/CNEMONSVariable.hpp +++ b/SU2_CFD/include/variables/CNEMONSVariable.hpp @@ -72,30 +72,12 @@ class CNEMONSVariable final : public CNEMOEulerVariable { * \param[in] val_nPrimVargrad - Number of primitive gradient variables. * \param[in] config - Definition of the particular problem. */ - CNEMONSVariable(su2double val_density, const su2double *val_massfrac, su2double *val_velocity, + CNEMONSVariable(su2double val_density, const su2double *val_massfrac, const su2double *val_velocity, su2double val_temperature, su2double val_temperature_ve, unsigned long npoint, unsigned long val_nDim, unsigned long val_nVar, unsigned long val_nPrimVar, - unsigned long val_nPrimVarGrad, CConfig *config, CNEMOGas *fluidmodel); + unsigned long val_nPrimVarGrad, const CConfig *config, CNEMOGas *fluidmodel); /*! - * \brief Constructor of the class. - * \param[in] val_solution - Pointer to the flow value (initialization value). - * \param[in] val_nDim - Number of dimensions of the problem. - * \param[in] val_nVar - Number of conserved variables. - * \param[in] val_nPrimVar - Number of primitive variables. - * \param[in] val_nPrimgVarGrad - Number of primitive gradient variables. - * \param[in] config - Definition of the particular problem. - */ - CNEMONSVariable(su2double *val_solution, unsigned long val_nDim, unsigned long val_nVar, - unsigned long val_nPrimVar, unsigned long val_nPrimVarGrad, unsigned long npoint, - CConfig *config); - - /*! - * \brief Destructor of the class. - */ - ~CNEMONSVariable() = default; - - /*! * \brief Get the primitive variables for all points. * \return Reference to primitives. */ @@ -106,11 +88,6 @@ class CNEMONSVariable final : public CNEMOEulerVariable { */ bool SetPrimVar(unsigned long iPoint, CFluidModel *FluidModel) final; - /*! - * \brief Set the vorticity value. - */ - bool SetVorticity(void); - /*! * \overload * \param[in] eddy_visc - Value of the eddy viscosity. diff --git a/SU2_CFD/include/variables/CNSVariable.hpp b/SU2_CFD/include/variables/CNSVariable.hpp index f1a5d93bc971..250910027d5f 100644 --- a/SU2_CFD/include/variables/CNSVariable.hpp +++ b/SU2_CFD/include/variables/CNSVariable.hpp @@ -56,12 +56,7 @@ class CNSVariable final : public CEulerVariable { * \param[in] config - Definition of the particular problem. */ CNSVariable(su2double density, const su2double *velocity, su2double energy, - unsigned long npoint, unsigned long ndim, unsigned long nvar, CConfig *config); - - /*! - * \brief Destructor of the class. - */ - ~CNSVariable() override = default; + unsigned long npoint, unsigned long ndim, unsigned long nvar, const CConfig *config); /*! * \brief Set the laminar viscosity. diff --git a/SU2_CFD/include/variables/CScalarVariable.hpp b/SU2_CFD/include/variables/CScalarVariable.hpp index cb4db4f1c33b..640373d60135 100644 --- a/SU2_CFD/include/variables/CScalarVariable.hpp +++ b/SU2_CFD/include/variables/CScalarVariable.hpp @@ -50,12 +50,7 @@ class CScalarVariable : public CVariable { * \param[in] nvar - Number of variables of the problem. * \param[in] config - Definition of the particular problem. */ - CScalarVariable(unsigned long npoint, unsigned long ndim, unsigned long nvar, CConfig* config); - - /*! - * \brief Destructor of the class. - */ - ~CScalarVariable() override = default; + CScalarVariable(unsigned long npoint, unsigned long ndim, unsigned long nvar, const CConfig* config); /*! * \brief Get the array of the reconstruction variables gradient at a node. diff --git a/SU2_CFD/include/variables/CVariable.hpp b/SU2_CFD/include/variables/CVariable.hpp index df3bba2712d1..7acc31f7fbec 100644 --- a/SU2_CFD/include/variables/CVariable.hpp +++ b/SU2_CFD/include/variables/CVariable.hpp @@ -142,7 +142,7 @@ class CVariable { * \param[in] nvar - Number of variables of the problem. * \param[in] config - Definition of the particular problem. */ - CVariable(unsigned long npoint, unsigned long nvar, CConfig *config); + CVariable(unsigned long npoint, unsigned long nvar, const CConfig *config); /*! * \overload @@ -152,7 +152,7 @@ class CVariable { * \param[in] config - Definition of the particular problem. * \param[in] adjoint - True if derived class is an adjoint variable. */ - CVariable(unsigned long npoint, unsigned long ndim, unsigned long nvar, CConfig *config, bool adjoint = false); + CVariable(unsigned long npoint, unsigned long ndim, unsigned long nvar, const CConfig *config, bool adjoint = false); /*! * \brief Destructor of the class. diff --git a/SU2_CFD/obj/Makefile.am b/SU2_CFD/obj/Makefile.am index 50950b471ad0..5980a61bf633 100644 --- a/SU2_CFD/obj/Makefile.am +++ b/SU2_CFD/obj/Makefile.am @@ -219,7 +219,8 @@ libSU2Core_sources = ../src/definition_structure.cpp \ ../src/variables/CAdjEulerVariable.cpp \ ../src/variables/CDiscAdjVariable.cpp \ ../src/variables/CIncNSVariable.cpp \ - ../src/variables/CEulerVariable.cpp + ../src/variables/CEulerVariable.cpp \ + ../src/variables/CFlowVariable.cpp su2_cfd_sources = \ ../src/SU2_CFD.cpp diff --git a/SU2_CFD/src/integration/CMultiGridIntegration.cpp b/SU2_CFD/src/integration/CMultiGridIntegration.cpp index 29794722b440..527c4beb6515 100644 --- a/SU2_CFD/src/integration/CMultiGridIntegration.cpp +++ b/SU2_CFD/src/integration/CMultiGridIntegration.cpp @@ -185,8 +185,7 @@ void CMultiGridIntegration::MultiGrid_Cycle(CGeometry ****geometry, solver_fine->Set_OldSolution(); - if (classical_rk4) - solver_fine->Set_NewSolution(); + if (classical_rk4) solver_fine->Set_NewSolution(); /*--- Compute time step, max eigenvalue, and integration scheme (steady and unsteady problems) ---*/ diff --git a/SU2_CFD/src/meson.build b/SU2_CFD/src/meson.build index 66f4ecbfde76..28a04702e764 100644 --- a/SU2_CFD/src/meson.build +++ b/SU2_CFD/src/meson.build @@ -72,6 +72,7 @@ su2_cfd_src += files(['variables/CIncNSVariable.cpp', 'variables/CVariable.cpp', 'variables/CNSVariable.cpp', 'variables/CAdjTurbVariable.cpp', + 'variables/CFlowVariable.cpp', 'variables/CIncEulerVariable.cpp', 'variables/CEulerVariable.cpp', 'variables/CNEMOEulerVariable.cpp', diff --git a/SU2_CFD/src/output/CNEMOCompOutput.cpp b/SU2_CFD/src/output/CNEMOCompOutput.cpp index f97f5cec3ed7..405f9c73a570 100644 --- a/SU2_CFD/src/output/CNEMOCompOutput.cpp +++ b/SU2_CFD/src/output/CNEMOCompOutput.cpp @@ -377,15 +377,17 @@ void CNEMOCompOutput::SetVolumeOutputFields(CConfig *config){ break; } - // Limiter values - AddVolumeOutput("LIMITER_DENSITY", "Limiter_Density", "LIMITER", "Limiter value of the density"); - AddVolumeOutput("LIMITER_MOMENTUM-X", "Limiter_Momentum_x", "LIMITER", "Limiter value of the x-momentum"); - AddVolumeOutput("LIMITER_MOMENTUM-Y", "Limiter_Momentum_y", "LIMITER", "Limiter value of the y-momentum"); - if (nDim == 3) - AddVolumeOutput("LIMITER_MOMENTUM-Z", "Limiter_Momentum_z", "LIMITER", "Limiter value of the z-momentum"); - AddVolumeOutput("LIMITER_ENERGY", "Limiter_Energy", "LIMITER", "Limiter value of the energy"); + if (config->GetKind_SlopeLimit_Flow() != NO_LIMITER && config->GetKind_SlopeLimit_Flow() != VAN_ALBADA_EDGE) { + // Limiter values + AddVolumeOutput("LIMITER_DENSITY", "Limiter_Density", "LIMITER", "Limiter value of the density"); + AddVolumeOutput("LIMITER_MOMENTUM-X", "Limiter_Momentum_x", "LIMITER", "Limiter value of the x-momentum"); + AddVolumeOutput("LIMITER_MOMENTUM-Y", "Limiter_Momentum_y", "LIMITER", "Limiter value of the y-momentum"); + if (nDim == 3) + AddVolumeOutput("LIMITER_MOMENTUM-Z", "Limiter_Momentum_z", "LIMITER", "Limiter value of the z-momentum"); + AddVolumeOutput("LIMITER_ENERGY", "Limiter_Energy", "LIMITER", "Limiter value of the energy"); + } - if (config->GetKind_SlopeLimit_Flow() != NO_LIMITER) { + if (config->GetKind_SlopeLimit_Turb() != NO_LIMITER) { switch(config->GetKind_Turb_Model()){ case SST: case SST_SUST: AddVolumeOutput("LIMITER_TKE", "Limiter_TKE", "LIMITER", "Limiter value of turb. kinetic energy"); @@ -519,17 +521,19 @@ void CNEMOCompOutput::LoadVolumeData(CConfig *config, CGeometry *geometry, CSolv break; } - SetVolumeOutputValue("LIMITER_DENSITY", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 0)); - SetVolumeOutputValue("LIMITER_MOMENTUM-X", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 1)); - SetVolumeOutputValue("LIMITER_MOMENTUM-Y", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 2)); - if (nDim == 3){ - SetVolumeOutputValue("LIMITER_MOMENTUM-Z", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 3)); - SetVolumeOutputValue("LIMITER_ENERGY", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 4)); - } else { - SetVolumeOutputValue("LIMITER_ENERGY", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 3)); + if (config->GetKind_SlopeLimit_Flow() != NO_LIMITER && config->GetKind_SlopeLimit_Flow() != VAN_ALBADA_EDGE) { + SetVolumeOutputValue("LIMITER_DENSITY", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 0)); + SetVolumeOutputValue("LIMITER_MOMENTUM-X", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 1)); + SetVolumeOutputValue("LIMITER_MOMENTUM-Y", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 2)); + if (nDim == 3){ + SetVolumeOutputValue("LIMITER_MOMENTUM-Z", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 3)); + SetVolumeOutputValue("LIMITER_ENERGY", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 4)); + } else { + SetVolumeOutputValue("LIMITER_ENERGY", iPoint, Node_Flow->GetLimiter_Primitive(iPoint, 3)); + } } - if (config->GetKind_SlopeLimit_Flow() != NO_LIMITER) { + if (config->GetKind_SlopeLimit_Turb() != NO_LIMITER) { switch(config->GetKind_Turb_Model()){ case SST: case SST_SUST: SetVolumeOutputValue("LIMITER_TKE", iPoint, Node_Turb->GetLimiter(iPoint, 0)); diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index aadcae81a3ea..2cd9ba5569ff 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -33,6 +33,7 @@ #include "../../include/fluid/CVanDerWaalsGas.hpp" #include "../../include/fluid/CPengRobinson.hpp" #include "../../include/numerics_simd/CNumericsSIMD.hpp" +#include "../../include/limiters/CLimiterDetails.hpp" CEulerSolver::CEulerSolver(CGeometry *geometry, CConfig *config, @@ -1782,8 +1783,8 @@ void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_contain if (van_albada) { su2double V_ij = V_j[iVar] - V_i[iVar]; - lim_i = V_ij*( 2.0*Project_Grad_i + V_ij) / (4*pow(Project_Grad_i, 2) + pow(V_ij, 2) + EPS); - lim_j = V_ij*(-2.0*Project_Grad_j + V_ij) / (4*pow(Project_Grad_j, 2) + pow(V_ij, 2) + EPS); + lim_i = LimiterHelpers<>::vanAlbadaFunction(Project_Grad_i, V_ij, EPS); + lim_j = LimiterHelpers<>::vanAlbadaFunction(-Project_Grad_j, V_ij, EPS); } else if (limiter) { lim_i = nodes->GetLimiter_Primitive(iPoint, iVar); diff --git a/SU2_CFD/src/solvers/CIncEulerSolver.cpp b/SU2_CFD/src/solvers/CIncEulerSolver.cpp index 2f210f27c3f0..52644844da0d 100644 --- a/SU2_CFD/src/solvers/CIncEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CIncEulerSolver.cpp @@ -31,6 +31,7 @@ #include "../../include/fluid/CIncIdealGas.hpp" #include "../../include/fluid/CIncIdealGasPolynomial.hpp" #include "../../include/variables/CIncNSVariable.hpp" +#include "../../include/limiters/CLimiterDetails.hpp" #include "../../../Common/include/toolboxes/geometry_toolbox.hpp" @@ -1184,8 +1185,8 @@ void CIncEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_cont if (van_albada) { su2double V_ij = V_j[iVar] - V_i[iVar]; - lim_i = V_ij*( 2.0*Project_Grad_i + V_ij) / (4*pow(Project_Grad_i, 2) + pow(V_ij, 2) + EPS); - lim_j = V_ij*(-2.0*Project_Grad_j + V_ij) / (4*pow(Project_Grad_j, 2) + pow(V_ij, 2) + EPS); + lim_i = LimiterHelpers<>::vanAlbadaFunction(Project_Grad_i, V_ij, EPS); + lim_j = LimiterHelpers<>::vanAlbadaFunction(-Project_Grad_j, V_ij, EPS); } else if (limiter) { lim_i = nodes->GetLimiter_Primitive(iPoint, iVar); diff --git a/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp b/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp index 2fda041c9bec..da99785d0547 100644 --- a/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp @@ -31,6 +31,7 @@ #include "../../../Common/include/toolboxes/printing_toolbox.hpp" #include "../../include/fluid/CMutationTCLib.hpp" #include "../../include/fluid/CSU2TCLib.hpp" +#include "../../include/limiters/CLimiterDetails.hpp" CNEMOEulerSolver::CNEMOEulerSolver(CGeometry *geometry, CConfig *config, unsigned short iMesh, const bool navier_stokes) : @@ -261,7 +262,7 @@ void CNEMOEulerSolver::CommonPreprocessing(CGeometry *geometry, CSolver **solver unsigned long tmp = ErrorCounter; SU2_MPI::Allreduce(&tmp, &ErrorCounter, 1, MPI_UNSIGNED_LONG, MPI_SUM, SU2_MPI::GetComm()); config->SetNonphysical_Points(ErrorCounter); - + if ((rank == MASTER_NODE) && (ErrorCounter != 0)) cout << "Warning. The initial solution contains "<< ErrorCounter << " points that are not physical." << endl; } @@ -547,18 +548,8 @@ void CNEMOEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_con auto Gradient_i = nodes->GetGradient_Reconstruction(iPoint); auto Gradient_j = nodes->GetGradient_Reconstruction(jPoint); - /*--- Set and extract limiters ---*/ - su2double *Limiter_i = nullptr, *Limiter_j = nullptr; - - if (limiter){ - Limiter_i = nodes->GetLimiter_Primitive(iPoint); - Limiter_j = nodes->GetLimiter_Primitive(jPoint); - } - for (iVar = 0; iVar < nPrimVarGrad; iVar++) { - su2double lim_i = 0.0; - su2double lim_j = 0.0; su2double Project_Grad_i = 0.0; su2double Project_Grad_j = 0.0; @@ -567,22 +558,22 @@ void CNEMOEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_con Project_Grad_j -= Vector_ij[iDim]*Gradient_j[iVar][iDim]; } - if (limiter) { - if (van_albada) { - su2double V_ij = V_j[iVar] - V_i[iVar]; - Limiter_i[iVar] = V_ij*( 2.0*Project_Grad_i + V_ij) / (4*pow(Project_Grad_i, 2) + pow(V_ij, 2) + EPS); - Limiter_j[iVar] = V_ij*(-2.0*Project_Grad_j + V_ij) / (4*pow(Project_Grad_j, 2) + pow(V_ij, 2) + EPS); - } - if (lim_i > Limiter_i[iVar] && Limiter_i[iVar] != 0) lim_i = Limiter_i[iVar]; - if (lim_j > Limiter_j[iVar] && Limiter_j[iVar] != 0) lim_j = Limiter_j[iVar]; - su2double lim_ij = min(lim_i, lim_j); - - Primitive_i[iVar] = V_i[iVar] + lim_ij*Project_Grad_i; - Primitive_j[iVar] = V_j[iVar] + lim_ij*Project_Grad_j; - } else { - Primitive_i[iVar] = V_i[iVar] + Project_Grad_i; - Primitive_j[iVar] = V_j[iVar] + Project_Grad_j; + su2double lim_i = 1.0; + su2double lim_j = 1.0; + + if (van_albada) { + su2double V_ij = V_j[iVar] - V_i[iVar]; + lim_i = LimiterHelpers<>::vanAlbadaFunction(Project_Grad_i, V_ij, EPS); + lim_j = LimiterHelpers<>::vanAlbadaFunction(-Project_Grad_j, V_ij, EPS); } + else if (limiter) { + lim_i = nodes->GetLimiter_Primitive(iPoint, iVar); + lim_j = nodes->GetLimiter_Primitive(jPoint, iVar); + } + su2double lim_ij = min(lim_i, lim_j); + + Primitive_i[iVar] = V_i[iVar] + lim_ij*Project_Grad_i; + Primitive_j[iVar] = V_j[iVar] + lim_ij*Project_Grad_j; } /*--- Check for non-physical solutions after reconstruction. If found, use the diff --git a/SU2_CFD/src/solvers/CTurbSASolver.cpp b/SU2_CFD/src/solvers/CTurbSASolver.cpp index 7beab143042e..c785b46861a0 100644 --- a/SU2_CFD/src/solvers/CTurbSASolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSASolver.cpp @@ -27,6 +27,7 @@ #include "../../include/solvers/CTurbSASolver.hpp" #include "../../include/variables/CTurbSAVariable.hpp" +#include "../../include/variables/CFlowVariable.hpp" #include "../../../Common/include/parallelization/omp_structure.hpp" #include "../../../Common/include/toolboxes/geometry_toolbox.hpp" @@ -226,11 +227,13 @@ void CTurbSASolver::Preprocessing(CGeometry *geometry, CSolver **solver_containe /*--- Set the vortex tilting coefficient at every node if required ---*/ if (kind_hybridRANSLES == SA_EDDES){ + auto* flowNodes = su2staticcast_p(solver_container[FLOW_SOL]->GetNodes()); + SU2_OMP_FOR_STAT(omp_chunk_size) for (unsigned long iPoint = 0; iPoint < nPoint; iPoint++){ - auto Vorticity = solver_container[FLOW_SOL]->GetNodes()->GetVorticity(iPoint); - auto PrimGrad_Flow = solver_container[FLOW_SOL]->GetNodes()->GetGradient_Primitive(iPoint); - auto Laminar_Viscosity = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); + auto Vorticity = flowNodes->GetVorticity(iPoint); + auto PrimGrad_Flow = flowNodes->GetGradient_Primitive(iPoint); + auto Laminar_Viscosity = flowNodes->GetLaminarViscosity(iPoint); nodes->SetVortex_Tilting(iPoint, PrimGrad_Flow, Vorticity, Laminar_Viscosity); } END_SU2_OMP_FOR @@ -309,10 +312,10 @@ void CTurbSASolver::Source_Residual(CGeometry *geometry, CSolver **solver_contai const bool transition = (config->GetKind_Trans_Model() == LM); const bool transition_BC = (config->GetKind_Trans_Model() == BC); - CVariable* flowNodes = solver_container[FLOW_SOL]->GetNodes(); + auto* flowNodes = su2staticcast_p(solver_container[FLOW_SOL]->GetNodes()); /*--- Pick one numerics object per thread. ---*/ - CNumerics* numerics = numerics_container[SOURCE_FIRST_TERM + omp_get_thread_num()*MAX_TERMS]; + auto* numerics = numerics_container[SOURCE_FIRST_TERM + omp_get_thread_num()*MAX_TERMS]; AD::StartNoSharedReading(); @@ -1347,42 +1350,33 @@ void CTurbSASolver::SetTurbVars_WF(CGeometry *geometry, CSolver **solver_contain void CTurbSASolver::SetDES_LengthScale(CSolver **solver, CGeometry *geometry, CConfig *config){ - unsigned short kindHybridRANSLES = config->GetKind_HybridRANSLES(); - unsigned long iPoint = 0, jPoint = 0; - unsigned short iDim = 0, jDim = 0, iNeigh = 0, nNeigh = 0; - - su2double constDES = config->GetConst_DES(); + const auto kindHybridRANSLES = config->GetKind_HybridRANSLES(); - su2double density = 0.0, laminarViscosity = 0.0, kinematicViscosity = 0.0, - eddyViscosity = 0.0, kinematicViscosityTurb = 0.0, wallDistance = 0.0, lengthScale = 0.0; + const su2double constDES = config->GetConst_DES(); - su2double maxDelta = 0.0, deltaAux = 0.0, distDES = 0.0, uijuij = 0.0, k2 = 0.0, r_d = 0.0, f_d = 0.0; - su2double deltaDDES = 0.0, omega = 0.0, ln_max = 0.0, ln[3] = {0.0}, aux_ln = 0.0, f_kh = 0.0; + const su2double fw_star = 0.424, cv1_3 = pow(7.1, 3), k2 = pow(0.41, 2); + const su2double cb1 = 0.1355, ct3 = 1.2, ct4 = 0.5; + const su2double sigma = 2./3., cb2 = 0.622, f_max = 1.0, f_min = 0.1, a1 = 0.15, a2 = 0.3; - su2double nu_hat, fw_star = 0.424, cv1_3 = pow(7.1, 3.0); k2 = pow(0.41, 2.0); - su2double cb1 = 0.1355, ct3 = 1.2, ct4 = 0.5; - su2double sigma = 2./3., cb2 = 0.622, f_max=1.0, f_min=0.1, a1=0.15, a2=0.3; - su2double cw1 = 0.0, Ji = 0.0, Ji_2 = 0.0, Ji_3 = 0.0, fv1 = 0.0, fv2 = 0.0, ft2 = 0.0, psi_2 = 0.0; - const su2double *coord_i = nullptr, *coord_j = nullptr, *vorticity = nullptr; - su2double delta[3] = {0.0}, ratioOmega[3] = {0.0}, vortexTiltingMeasure = 0.0; + auto* flowNodes = su2staticcast_p(solver[FLOW_SOL]->GetNodes()); SU2_OMP_FOR_DYN(omp_chunk_size) - for (iPoint = 0; iPoint < nPointDomain; iPoint++){ - - coord_i = geometry->nodes->GetCoord(iPoint); - nNeigh = geometry->nodes->GetnPoint(iPoint); - wallDistance = geometry->nodes->GetWall_Distance(iPoint); - const auto primVarGrad = solver[FLOW_SOL]->GetNodes()->GetGradient_Primitive(iPoint); - vorticity = solver[FLOW_SOL]->GetNodes()->GetVorticity(iPoint); - density = solver[FLOW_SOL]->GetNodes()->GetDensity(iPoint); - laminarViscosity = solver[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); - eddyViscosity = solver[TURB_SOL]->GetNodes()->GetmuT(iPoint); - kinematicViscosity = laminarViscosity/density; - kinematicViscosityTurb = eddyViscosity/density; - - uijuij = 0.0; - for(iDim = 0; iDim < nDim; iDim++){ - for(jDim = 0; jDim < nDim; jDim++){ + for (auto iPoint = 0ul; iPoint < nPointDomain; iPoint++){ + + const auto coord_i = geometry->nodes->GetCoord(iPoint); + const auto nNeigh = geometry->nodes->GetnPoint(iPoint); + const auto wallDistance = geometry->nodes->GetWall_Distance(iPoint); + const auto primVarGrad = flowNodes->GetGradient_Primitive(iPoint); + const auto vorticity = flowNodes->GetVorticity(iPoint); + const auto density = flowNodes->GetDensity(iPoint); + const auto laminarViscosity = flowNodes->GetLaminarViscosity(iPoint); + const auto eddyViscosity = nodes->GetmuT(iPoint); + const su2double kinematicViscosity = laminarViscosity/density; + const su2double kinematicViscosityTurb = eddyViscosity/density; + + su2double uijuij = 0.0; + for(auto iDim = 0u; iDim < nDim; iDim++){ + for(auto jDim = 0u; jDim < nDim; jDim++){ uijuij += primVarGrad[1+iDim][jDim]*primVarGrad[1+iDim][jDim]; } } @@ -1391,137 +1385,141 @@ void CTurbSASolver::SetDES_LengthScale(CSolver **solver, CGeometry *geometry, CC /*--- Low Reynolds number correction term ---*/ - nu_hat = nodes->GetSolution(iPoint,0); - Ji = nu_hat/kinematicViscosity; - Ji_2 = Ji * Ji; - Ji_3 = Ji*Ji*Ji; - fv1 = Ji_3/(Ji_3+cv1_3); - fv2 = 1.0 - Ji/(1.0+Ji*fv1); - ft2 = ct3*exp(-ct4*Ji_2); - cw1 = cb1/k2+(1.0+cb2)/sigma; + const su2double nu_hat = nodes->GetSolution(iPoint,0); + const su2double Ji = nu_hat/kinematicViscosity; + const su2double Ji_2 = Ji * Ji; + const su2double Ji_3 = Ji*Ji*Ji; + const su2double fv1 = Ji_3/(Ji_3+cv1_3); + const su2double fv2 = 1.0 - Ji/(1.0+Ji*fv1); + const su2double ft2 = ct3*exp(-ct4*Ji_2); + const su2double cw1 = cb1/k2+(1.0+cb2)/sigma; - psi_2 = (1.0 - (cb1/(cw1*k2*fw_star))*(ft2 + (1.0 - ft2)*fv2))/(fv1 * max(1.0e-10,1.0-ft2)); + su2double psi_2 = (1.0 - (cb1/(cw1*k2*fw_star))*(ft2 + (1.0 - ft2)*fv2))/(fv1 * max(1.0e-10,1.0-ft2)); psi_2 = min(100.0,psi_2); + su2double lengthScale = 0.0; + switch(kindHybridRANSLES){ - case SA_DES: + case SA_DES: { /*--- Original Detached Eddy Simulation (DES97) Spalart 1997 ---*/ - maxDelta = geometry->nodes->GetMaxLength(iPoint); - distDES = constDES * maxDelta; + const su2double maxDelta = geometry->nodes->GetMaxLength(iPoint); + const su2double distDES = constDES * maxDelta; lengthScale = min(distDES,wallDistance); break; - - case SA_DDES: + } + case SA_DDES: { /*--- A New Version of Detached-eddy Simulation, Resistant to Ambiguous Grid Densities. Spalart et al. Theoretical and Computational Fluid Dynamics - 2006 ---*/ - maxDelta = geometry->nodes->GetMaxLength(iPoint); + const su2double maxDelta = geometry->nodes->GetMaxLength(iPoint); - r_d = (kinematicViscosityTurb+kinematicViscosity)/(uijuij*k2*pow(wallDistance, 2.0)); - f_d = 1.0-tanh(pow(8.0*r_d,3.0)); + const su2double r_d = (kinematicViscosityTurb+kinematicViscosity)/(uijuij*k2*pow(wallDistance, 2.0)); + const su2double f_d = 1.0-tanh(pow(8.0*r_d,3.0)); - distDES = constDES * maxDelta; + const su2double distDES = constDES * maxDelta; lengthScale = wallDistance-f_d*max(0.0,(wallDistance-distDES)); break; - case SA_ZDES: + } + case SA_ZDES: { /*--- Recent improvements in the Zonal Detached Eddy Simulation (ZDES) formulation. Deck Theoretical and Computational Fluid Dynamics - 2012 ---*/ - for (iNeigh = 0; iNeigh < nNeigh; iNeigh++){ - jPoint = geometry->nodes->GetPoint(iPoint, iNeigh); - coord_j = geometry->nodes->GetCoord(jPoint); - for ( iDim = 0; iDim < nDim; iDim++){ - deltaAux = abs(coord_j[iDim] - coord_i[iDim]); - delta[iDim] = max(delta[iDim], deltaAux); - } - deltaDDES = geometry->nodes->GetMaxLength(iPoint); + const su2double deltaDDES = geometry->nodes->GetMaxLength(iPoint); + + su2double delta[MAXNDIM] = {}, ratioOmega[MAXNDIM] = {}; + + for (const auto jPoint : geometry->nodes->GetPoints(iPoint)) { + const auto coord_j = geometry->nodes->GetCoord(jPoint); + for (auto iDim = 0u; iDim < nDim; iDim++){ + const su2double deltaAux = abs(coord_j[iDim] - coord_i[iDim]); + delta[iDim] = max(delta[iDim], deltaAux); + } } - omega = sqrt(vorticity[0]*vorticity[0] + - vorticity[1]*vorticity[1] + - vorticity[2]*vorticity[2]); + const su2double omega = GeometryToolbox::Norm(3, vorticity); - for (iDim = 0; iDim < 3; iDim++){ + for (auto iDim = 0u; iDim < 3; iDim++){ ratioOmega[iDim] = vorticity[iDim]/omega; } - maxDelta = sqrt(pow(ratioOmega[0],2.0)*delta[1]*delta[2] + - pow(ratioOmega[1],2.0)*delta[0]*delta[2] + - pow(ratioOmega[2],2.0)*delta[0]*delta[1]); + su2double maxDelta = sqrt(pow(ratioOmega[0], 2)*delta[1]*delta[2] + + pow(ratioOmega[1], 2)*delta[0]*delta[2] + + pow(ratioOmega[2], 2)*delta[0]*delta[1]); - r_d = (kinematicViscosityTurb+kinematicViscosity)/(uijuij*k2*pow(wallDistance, 2.0)); - f_d = 1.0-tanh(pow(8.0*r_d,3.0)); + const su2double r_d = (kinematicViscosityTurb+kinematicViscosity)/(uijuij*k2*pow(wallDistance, 2.0)); + const su2double f_d = 1.0-tanh(pow(8.0*r_d,3.0)); if (f_d < 0.99){ maxDelta = deltaDDES; } - distDES = constDES * maxDelta; + const su2double distDES = constDES * maxDelta; lengthScale = wallDistance-f_d*max(0.0,(wallDistance-distDES)); break; - - case SA_EDDES: - + } + case SA_EDDES: { /*--- An Enhanced Version of DES with Rapid Transition from RANS to LES in Separated Flows. Shur et al. Flow Turbulence Combust - 2015 ---*/ - vortexTiltingMeasure = nodes->GetVortex_Tilting(iPoint); + su2double vortexTiltingMeasure = nodes->GetVortex_Tilting(iPoint); - omega = sqrt(vorticity[0]*vorticity[0] + - vorticity[1]*vorticity[1] + - vorticity[2]*vorticity[2]); + const su2double omega = GeometryToolbox::Norm(3, vorticity); - for (iDim = 0; iDim < 3; iDim++){ + su2double ratioOmega[MAXNDIM] = {}; + + for (auto iDim = 0; iDim < 3; iDim++){ ratioOmega[iDim] = vorticity[iDim]/omega; } - ln_max = 0.0; - deltaDDES = 0.0; - for (iNeigh = 0;iNeigh < nNeigh; iNeigh++){ - jPoint = geometry->nodes->GetPoint(iPoint, iNeigh); - coord_j = geometry->nodes->GetCoord(jPoint); - for (iDim = 0; iDim < nDim; iDim++){ + const su2double deltaDDES = geometry->nodes->GetMaxLength(iPoint); + + su2double ln_max = 0.0; + for (const auto jPoint : geometry->nodes->GetPoints(iPoint)) { + const auto coord_j = geometry->nodes->GetCoord(jPoint); + su2double delta[MAXNDIM] = {}; + for (auto iDim = 0u; iDim < nDim; iDim++){ delta[iDim] = fabs(coord_j[iDim] - coord_i[iDim]); } - deltaDDES = geometry->nodes->GetMaxLength(iPoint); + su2double ln[3]; ln[0] = delta[1]*ratioOmega[2] - delta[2]*ratioOmega[1]; ln[1] = delta[2]*ratioOmega[0] - delta[0]*ratioOmega[2]; ln[2] = delta[0]*ratioOmega[1] - delta[1]*ratioOmega[0]; - aux_ln = sqrt(ln[0]*ln[0] + ln[1]*ln[1] + ln[2]*ln[2]); - ln_max = max(ln_max,aux_ln); + const su2double aux_ln = sqrt(ln[0]*ln[0] + ln[1]*ln[1] + ln[2]*ln[2]); + ln_max = max(ln_max, aux_ln); vortexTiltingMeasure += nodes->GetVortex_Tilting(jPoint); } + vortexTiltingMeasure /= (nNeigh + 1); - vortexTiltingMeasure = (vortexTiltingMeasure/fabs(nNeigh + 1.0)); + const su2double f_kh = max(f_min, + min(f_max, + f_min + ((f_max - f_min)/(a2 - a1)) * (vortexTiltingMeasure - a1))); - f_kh = max(f_min, min(f_max, f_min + ((f_max - f_min)/(a2 - a1)) * (vortexTiltingMeasure - a1))); + const su2double r_d = (kinematicViscosityTurb+kinematicViscosity)/(uijuij*k2*pow(wallDistance, 2.0)); + const su2double f_d = 1.0-tanh(pow(8.0*r_d,3.0)); - r_d = (kinematicViscosityTurb+kinematicViscosity)/(uijuij*k2*pow(wallDistance, 2.0)); - f_d = 1.0-tanh(pow(8.0*r_d,3.0)); - - maxDelta = (ln_max/sqrt(3.0)) * f_kh; + su2double maxDelta = (ln_max/sqrt(3.0)) * f_kh; if (f_d < 0.999){ maxDelta = deltaDDES; } - distDES = constDES * maxDelta; - lengthScale=wallDistance-f_d*max(0.0,(wallDistance-distDES)); + const su2double distDES = constDES * maxDelta; + lengthScale = wallDistance-f_d*max(0.0,(wallDistance-distDES)); break; - + } } nodes->SetDES_LengthScale(iPoint, lengthScale); diff --git a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp index 3245a65e24ce..97e18236fdbe 100644 --- a/SU2_CFD/src/solvers/CTurbSSTSolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSSTSolver.cpp @@ -27,6 +27,7 @@ #include "../../include/solvers/CTurbSSTSolver.hpp" #include "../../include/variables/CTurbSSTVariable.hpp" +#include "../../include/variables/CFlowVariable.hpp" #include "../../../Common/include/parallelization/omp_structure.hpp" #include "../../../Common/include/toolboxes/geometry_toolbox.hpp" @@ -241,20 +242,19 @@ void CTurbSSTSolver::Postprocessing(CGeometry *geometry, CSolver **solver_contai AD::StartNoSharedReading(); + auto* flowNodes = su2staticcast_p(solver_container[FLOW_SOL]->GetNodes()); + SU2_OMP_FOR_STAT(omp_chunk_size) for (unsigned long iPoint = 0; iPoint < nPoint; iPoint ++) { /*--- Compute blending functions and cross diffusion ---*/ - su2double rho = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); - su2double mu = solver_container[FLOW_SOL]->GetNodes()->GetLaminarViscosity(iPoint); + su2double rho = flowNodes->GetDensity(iPoint); + su2double mu = flowNodes->GetLaminarViscosity(iPoint); su2double dist = geometry->nodes->GetWall_Distance(iPoint); - const su2double *Vorticity = solver_container[FLOW_SOL]->GetNodes()->GetVorticity(iPoint); - su2double VorticityMag = sqrt(Vorticity[0]*Vorticity[0] + - Vorticity[1]*Vorticity[1] + - Vorticity[2]*Vorticity[2]); + su2double VorticityMag = GeometryToolbox::Norm(3, flowNodes->GetVorticity(iPoint)); nodes->SetBlendingFunc(iPoint, mu, dist, rho); @@ -296,10 +296,10 @@ void CTurbSSTSolver::Source_Residual(CGeometry *geometry, CSolver **solver_conta const bool implicit = (config->GetKind_TimeIntScheme() == EULER_IMPLICIT); - CVariable* flowNodes = solver_container[FLOW_SOL]->GetNodes(); + auto* flowNodes = su2staticcast_p(solver_container[FLOW_SOL]->GetNodes()); /*--- Pick one numerics object per thread. ---*/ - CNumerics* numerics = numerics_container[SOURCE_FIRST_TERM + omp_get_thread_num()*MAX_TERMS]; + auto* numerics = numerics_container[SOURCE_FIRST_TERM + omp_get_thread_num()*MAX_TERMS]; /*--- Loop over all points. ---*/ diff --git a/SU2_CFD/src/variables/CEulerVariable.cpp b/SU2_CFD/src/variables/CEulerVariable.cpp index 899a762cd593..1a493c3f74ba 100644 --- a/SU2_CFD/src/variables/CEulerVariable.cpp +++ b/SU2_CFD/src/variables/CEulerVariable.cpp @@ -29,50 +29,16 @@ #include "../../include/fluid/CFluidModel.hpp" CEulerVariable::CEulerVariable(su2double density, const su2double *velocity, su2double energy, unsigned long npoint, - unsigned long ndim, unsigned long nvar, CConfig *config) : CVariable(npoint, ndim, nvar, config), - Gradient_Reconstruction(config->GetReconstructionGradientRequired() ? Gradient_Aux : Gradient_Primitive) { + unsigned long ndim, unsigned long nvar, const CConfig *config) + : CFlowVariable(npoint, ndim, nvar, ndim + 9, ndim + 4, config) { const bool dual_time = (config->GetTime_Marching() == TIME_MARCHING::DT_STEPPING_1ST) || (config->GetTime_Marching() == TIME_MARCHING::DT_STEPPING_2ND); - const bool viscous = config->GetViscous(); - const bool windgust = config->GetWind_Gust(); const bool classical_rk4 = (config->GetKind_TimeIntScheme_Flow() == CLASSICAL_RK4_EXPLICIT); - /*--- Allocate and initialize the primitive variables and gradients ---*/ - - nPrimVar = nDim+9; - nPrimVarGrad = nDim+4; - nSecondaryVar = viscous? 8 : 2; + nSecondaryVar = config->GetViscous() ? 8 : 2, nSecondaryVarGrad = 2; - /*--- Allocate residual structures ---*/ - - Res_TruncError.resize(nPoint,nVar) = su2double(0.0); - - /*--- Only for residual smoothing (multigrid) ---*/ - - for (unsigned long iMesh = 0; iMesh <= config->GetnMGLevels(); iMesh++) { - if (config->GetMG_CorrecSmooth(iMesh) > 0) { - Residual_Sum.resize(nPoint,nVar); - Residual_Old.resize(nPoint,nVar); - break; - } - } - - /*--- Allocate undivided laplacian (centered) ---*/ - - if (config->GetKind_ConvNumScheme_Flow() == SPACE_CENTERED) - Undivided_Laplacian.resize(nPoint,nVar); - - /*--- Allocate the slope limiter (MUSCL upwind) ---*/ - - if (config->GetKind_SlopeLimit_Flow() != NO_LIMITER && - config->GetKind_SlopeLimit_Flow() != VAN_ALBADA_EDGE) { - Limiter_Primitive.resize(nPoint,nPrimVarGrad) = su2double(0.0); - Solution_Max.resize(nPoint,nPrimVarGrad) = su2double(0.0); - Solution_Min.resize(nPoint,nPrimVarGrad) = su2double(0.0); - } - /*--- Solution initialization ---*/ su2double val_solution[5] = {su2double(1.0), velocity[0], velocity[1], energy, energy}; @@ -84,8 +50,6 @@ CEulerVariable::CEulerVariable(su2double density, const su2double *velocity, su2 Solution_Old = Solution; - /*--- New solution initialization for Classical RK4 ---*/ - if (classical_rk4) Solution_New = Solution; /*--- Allocate and initializate solution for dual time strategy ---*/ @@ -95,61 +59,18 @@ CEulerVariable::CEulerVariable(su2double density, const su2double *velocity, su2 Solution_time_n1 = Solution; } - /*--- Allocate space for the harmonic balance source terms ---*/ - - if (config->GetTime_Marching() == TIME_MARCHING::HARMONIC_BALANCE) - HB_Source.resize(nPoint,nVar) = su2double(0.0); - - /*--- Allocate vector for wind gust and wind gust derivative field ---*/ - - if (windgust) { - WindGust.resize(nPoint,nDim); - WindGustDer.resize(nPoint,nDim+1); - } - - /*--- Compressible flow, primitive variables (T, vx, vy, vz, P, rho, h, c, mu, mut, k, Cp) ---*/ - - Primitive.resize(nPoint,nPrimVar) = su2double(0.0); Secondary.resize(nPoint,nSecondaryVar) = su2double(0.0); - /*--- Compressible flow, gradients primitive variables (T, vx, vy, vz, P, rho, h) ---*/ - - if (config->GetMUSCL_Flow() || viscous) { - Gradient_Primitive.resize(nPoint,nPrimVarGrad,nDim,0.0); - } - - if (config->GetReconstructionGradientRequired() && - config->GetKind_ConvNumScheme_Flow() != SPACE_CENTERED) { - Gradient_Aux.resize(nPoint,nPrimVarGrad,nDim,0.0); - } - if (config->GetAxisymmetric()){ nAuxVar = 3; Grad_AuxVar.resize(nPoint,nAuxVar,nDim,0.0); AuxVar.resize(nPoint,nAuxVar) = su2double(0.0); } - if (config->GetLeastSquaresRequired()) { - Rmatrix.resize(nPoint,nDim,nDim,0.0); + if (config->GetWind_Gust()) { + WindGust.resize(nPoint,nDim); + WindGustDer.resize(nPoint,nDim+1); } - - if (config->GetMultizone_Problem()) - Set_BGSSolution_k(); - - Velocity2.resize(nPoint) = su2double(0.0); - Max_Lambda_Inv.resize(nPoint) = su2double(0.0); - Delta_Time.resize(nPoint) = su2double(0.0); - Lambda.resize(nPoint) = su2double(0.0); - Sensor.resize(nPoint) = su2double(0.0); - - /* Under-relaxation parameter. */ - UnderRelaxation.resize(nPoint) = su2double(1.0); - LocalCFL.resize(nPoint) = su2double(0.0); - - /* Non-physical point (first-order) initialization. */ - Non_Physical.resize(nPoint) = false; - Non_Physical_Counter.resize(nPoint) = 0; - } bool CEulerVariable::SetPrimVar(unsigned long iPoint, CFluidModel *FluidModel) { @@ -208,5 +129,3 @@ void CEulerVariable::SetSecondaryVar(unsigned long iPoint, CFluidModel *FluidMod SetdPde_rho(iPoint, FluidModel->GetdPde_rho()); } - -void CEulerVariable::SetSolution_New() { Solution_New = Solution; } diff --git a/SU2_CFD/src/variables/CFlowVariable.cpp b/SU2_CFD/src/variables/CFlowVariable.cpp new file mode 100644 index 000000000000..16b5d55eb545 --- /dev/null +++ b/SU2_CFD/src/variables/CFlowVariable.cpp @@ -0,0 +1,105 @@ +/*! + * \file CFlowVariable.cpp + * \brief Definition of common solution fields for flow solvers. + * \version 7.2.0 "Blackbird" + * + * SU2 Project Website: https://su2code.github.io + * + * The SU2 Project is maintained by the SU2 Foundation + * (http://su2foundation.org) + * + * Copyright 2012-2021, 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 "../../include/variables/CFlowVariable.hpp" +#include "../../../Common/include/parallelization/omp_structure.hpp" + +CFlowVariable::CFlowVariable(unsigned long npoint, unsigned long ndim, unsigned long nvar, unsigned long nprimvar, + unsigned long nprimvargrad, const CConfig* config) + : CVariable(npoint, ndim, nvar, config), + Gradient_Reconstruction(config->GetReconstructionGradientRequired() ? Gradient_Aux : Gradient_Primitive) { + nPrimVar = nprimvar; + nPrimVarGrad = nprimvargrad; + + /*--- Allocate residual structures for multigrid. ---*/ + + Res_TruncError.resize(nPoint, nVar) = su2double(0.0); + + for (unsigned long iMesh = 0; iMesh <= config->GetnMGLevels(); iMesh++) { + if (config->GetMG_CorrecSmooth(iMesh) > 0) { + Residual_Sum.resize(nPoint, nVar); + Residual_Old.resize(nPoint, nVar); + break; + } + } + + /*--- Primitive variables and gradients (see derived classes for what is stored each column) ---*/ + + Primitive.resize(nPoint, nPrimVar) = su2double(0.0); + + if (config->GetMUSCL_Flow() || config->GetViscous()) { + Gradient_Primitive.resize(nPoint, nPrimVarGrad, nDim, 0.0); + } + + if (config->GetReconstructionGradientRequired() && config->GetKind_ConvNumScheme_Flow() != SPACE_CENTERED) { + Gradient_Aux.resize(nPoint, nPrimVarGrad, nDim, 0.0); + } + + if (config->GetLeastSquaresRequired()) { + Rmatrix.resize(nPoint, nDim, nDim, 0.0); + } + + /*--- Allocate undivided laplacian (centered) ---*/ + + if (config->GetKind_ConvNumScheme_Flow() == SPACE_CENTERED) { + Undivided_Laplacian.resize(nPoint, nVar); + } + + /*--- Allocate the slope limiter (MUSCL upwind) ---*/ + + if (config->GetKind_SlopeLimit_Flow() != NO_LIMITER && config->GetKind_SlopeLimit_Flow() != VAN_ALBADA_EDGE) { + Limiter_Primitive.resize(nPoint, nPrimVarGrad) = su2double(0.0); + Solution_Max.resize(nPoint, nPrimVarGrad) = su2double(0.0); + Solution_Min.resize(nPoint, nPrimVarGrad) = su2double(0.0); + } + + Velocity2.resize(nPoint) = su2double(0.0); + Max_Lambda_Inv.resize(nPoint) = su2double(0.0); + Delta_Time.resize(nPoint) = su2double(0.0); + Lambda.resize(nPoint) = su2double(0.0); + Sensor.resize(nPoint) = su2double(0.0); + + /* Under-relaxation parameter. */ + UnderRelaxation.resize(nPoint) = su2double(1.0); + LocalCFL.resize(nPoint) = su2double(0.0); + + /* Non-physical point (first-order) initialization. */ + Non_Physical.resize(nPoint) = false; + Non_Physical_Counter.resize(nPoint) = 0; + + if (config->GetMultizone_Problem()) { + Set_BGSSolution_k(); + } + + if (config->GetTime_Marching() == TIME_MARCHING::HARMONIC_BALANCE) { + HB_Source.resize(nPoint, nVar) = su2double(0.0); + } +} + +void CFlowVariable::SetSolution_New() { + assert(Solution_New.size() == Solution.size()); + parallelCopy(Solution.size(), Solution.data(), Solution_New.data()); +} diff --git a/SU2_CFD/src/variables/CIncEulerVariable.cpp b/SU2_CFD/src/variables/CIncEulerVariable.cpp index a181a2b918c8..dab0825899a2 100644 --- a/SU2_CFD/src/variables/CIncEulerVariable.cpp +++ b/SU2_CFD/src/variables/CIncEulerVariable.cpp @@ -28,46 +28,13 @@ #include "../../include/variables/CIncEulerVariable.hpp" #include "../../include/fluid/CFluidModel.hpp" -CIncEulerVariable::CIncEulerVariable(su2double pressure, const su2double *velocity, su2double temperature, unsigned long npoint, - unsigned long ndim, unsigned long nvar, CConfig *config) : CVariable(npoint, ndim, nvar, config), - Gradient_Reconstruction(config->GetReconstructionGradientRequired() ? Gradient_Aux : Gradient_Primitive) { +CIncEulerVariable::CIncEulerVariable(su2double pressure, const su2double *velocity, su2double temperature, + unsigned long npoint, unsigned long ndim, unsigned long nvar, const CConfig *config) + : CFlowVariable(npoint, ndim, nvar, ndim + 9, ndim + 4, config) { const bool dual_time = (config->GetTime_Marching() == TIME_MARCHING::DT_STEPPING_1ST) || (config->GetTime_Marching() == TIME_MARCHING::DT_STEPPING_2ND); - const bool viscous = config->GetViscous(); - - /*--- Allocate and initialize the primitive variables and gradients. - Make sure to align the sizes with the constructor of CIncEulerSolver ---*/ - - nPrimVar = nDim+9; nPrimVarGrad = nDim+4; - - /*--- Allocate residual structures ---*/ - - Res_TruncError.resize(nPoint,nVar) = su2double(0.0); - - /*--- Only for residual smoothing (multigrid) ---*/ - - for (unsigned long iMesh = 0; iMesh <= config->GetnMGLevels(); iMesh++) { - if (config->GetMG_CorrecSmooth(iMesh) > 0) { - Residual_Sum.resize(nPoint,nVar); - Residual_Old.resize(nPoint,nVar); - break; - } - } - - /*--- Allocate undivided laplacian (centered) ---*/ - - if (config->GetKind_ConvNumScheme_Flow() == SPACE_CENTERED) - Undivided_Laplacian.resize(nPoint,nVar); - - /*--- Allocate the slope limiter (MUSCL upwind) ---*/ - - if (config->GetKind_SlopeLimit_Flow() != NO_LIMITER && - config->GetKind_SlopeLimit_Flow() != VAN_ALBADA_EDGE) { - Limiter_Primitive.resize(nPoint,nPrimVarGrad) = su2double(0.0); - Solution_Max.resize(nPoint,nPrimVarGrad) = su2double(0.0); - Solution_Min.resize(nPoint,nPrimVarGrad) = su2double(0.0); - } + const bool classical_rk4 = (config->GetKind_TimeIntScheme_Flow() == CLASSICAL_RK4_EXPLICIT); /*--- Solution initialization ---*/ @@ -80,6 +47,8 @@ CIncEulerVariable::CIncEulerVariable(su2double pressure, const su2double *veloci Solution_Old = Solution; + if (classical_rk4) Solution_New = Solution; + /*--- Allocate and initialize solution for dual time strategy ---*/ if (dual_time) { @@ -87,48 +56,11 @@ CIncEulerVariable::CIncEulerVariable(su2double pressure, const su2double *veloci Solution_time_n1 = Solution; } - /*--- Incompressible flow, primitive variables nDim+9, (P, vx, vy, vz, T, rho, beta, lamMu, EddyMu, Kt_eff, Cp, Cv) ---*/ - - Primitive.resize(nPoint,nPrimVar) = su2double(0.0); - - /*--- Incompressible flow, gradients primitive variables nDim+4, (P, vx, vy, vz, T, rho, beta) ---*/ - - if (config->GetMUSCL_Flow() || viscous) { - Gradient_Primitive.resize(nPoint,nPrimVarGrad,nDim,0.0); - } - - if (config->GetReconstructionGradientRequired() && - config->GetKind_ConvNumScheme_Flow() != SPACE_CENTERED) { - Gradient_Aux.resize(nPoint,nPrimVarGrad,nDim,0.0); - } - - if (config->GetLeastSquaresRequired()) { - Rmatrix.resize(nPoint,nDim,nDim,0.0); - } - - if (config->GetMultizone_Problem()) - Set_BGSSolution_k(); - - Velocity2.resize(nPoint) = su2double(0.0); - Max_Lambda_Inv.resize(nPoint) = su2double(0.0); - Delta_Time.resize(nPoint) = su2double(0.0); - Lambda.resize(nPoint) = su2double(0.0); - Sensor.resize(nPoint) = su2double(0.0); - if (config->GetKind_Streamwise_Periodic() != ENUM_STREAMWISE_PERIODIC::NONE) { Streamwise_Periodic_RecoveredPressure.resize(nPoint) = su2double(0.0); if (config->GetStreamwise_Periodic_Temperature()) Streamwise_Periodic_RecoveredTemperature.resize(nPoint) = su2double(0.0); } - - /* Under-relaxation parameter. */ - UnderRelaxation.resize(nPoint) = su2double(1.0); - LocalCFL.resize(nPoint) = su2double(0.0); - - /* Non-physical point (first-order) initialization. */ - Non_Physical.resize(nPoint) = false; - Non_Physical_Counter.resize(nPoint) = 0; - } bool CIncEulerVariable::SetPrimVar(unsigned long iPoint, CFluidModel *FluidModel) { diff --git a/SU2_CFD/src/variables/CIncNSVariable.cpp b/SU2_CFD/src/variables/CIncNSVariable.cpp index 2f0c213820c1..60352af4b0ee 100644 --- a/SU2_CFD/src/variables/CIncNSVariable.cpp +++ b/SU2_CFD/src/variables/CIncNSVariable.cpp @@ -29,7 +29,7 @@ #include "../../include/fluid/CFluidModel.hpp" CIncNSVariable::CIncNSVariable(su2double pressure, const su2double *velocity, su2double temperature, - unsigned long npoint, unsigned long ndim, unsigned long nvar, CConfig *config) : + unsigned long npoint, unsigned long ndim, unsigned long nvar, const CConfig *config) : CIncEulerVariable(pressure, velocity, temperature, npoint, ndim, nvar, config) { Vorticity.resize(nPoint,3); @@ -37,14 +37,11 @@ CIncNSVariable::CIncNSVariable(su2double pressure, const su2double *velocity, su DES_LengthScale.resize(nPoint) = su2double(0.0); Max_Lambda_Visc.resize(nPoint); - if (config->GetAxisymmetric()) { - nAuxVar = 1; - AuxVar.resize(nPoint,nAuxVar) = su2double(0.0); - Grad_AuxVar.resize(nPoint,nAuxVar,nDim); - } - - /*--- Allocate memory for the AuxVar+gradient of eddy viscosity mu_t ---*/ - if (config->GetStreamwise_Periodic_Temperature() && (config->GetKind_Turb_Model() != NONE)) { + /*--- Allocate memory for the AuxVar and its gradient. See e.g. CIncEulerSolver::Source_Residual: + * Axisymmetric: total-viscosity * y-vel / y-coord + * Streamwise Periodic: eddy viscosity (mu_t) ---*/ + if (config->GetAxisymmetric() || + (config->GetStreamwise_Periodic_Temperature() && (config->GetKind_Turb_Model() != NONE))) { nAuxVar = 1; AuxVar.resize(nPoint,nAuxVar) = su2double(0.0); Grad_AuxVar.resize(nPoint,nAuxVar,nDim); diff --git a/SU2_CFD/src/variables/CNEMOEulerVariable.cpp b/SU2_CFD/src/variables/CNEMOEulerVariable.cpp index 562d1c1a70ba..d8589b1f79fd 100644 --- a/SU2_CFD/src/variables/CNEMOEulerVariable.cpp +++ b/SU2_CFD/src/variables/CNEMOEulerVariable.cpp @@ -30,7 +30,7 @@ CNEMOEulerVariable::CNEMOEulerVariable(su2double val_pressure, const su2double *val_massfrac, - su2double *val_mach, + const su2double *val_mach, su2double val_temperature, su2double val_temperature_ve, unsigned long npoint, @@ -38,20 +38,14 @@ CNEMOEulerVariable::CNEMOEulerVariable(su2double val_pressure, unsigned long nvar, unsigned long nvarprim, unsigned long nvarprimgrad, - CConfig *config, - CNEMOGas *fluidmodel) : CVariable(npoint, - ndim, - nvar, - config ), - Gradient_Reconstruction(config->GetReconstructionGradientRequired() ? Gradient_Aux : Gradient_Primitive), - implicit(config->GetKind_TimeIntScheme_Flow() == EULER_IMPLICIT) { + const CConfig *config, + CNEMOGas *fluidmodel) + : CFlowVariable(npoint, ndim, nvar, nvarprim, nvarprimgrad, config), + implicit(config->GetKind_TimeIntScheme_Flow() == EULER_IMPLICIT) { unsigned short iDim, iSpecies; /*--- Setting variable amounts ---*/ - nDim = ndim; - nPrimVar = nvarprim; - nPrimVarGrad = nvarprimgrad; nSpecies = config->GetnSpecies(); RHOS_INDEX = 0; @@ -73,9 +67,6 @@ CNEMOEulerVariable::CNEMOEulerVariable(su2double val_pressure, Tve_Freestream = config->GetTemperature_ve_FreeStream(); } - /*--- Allocate & initialize residual vectors ---*/ - Res_TruncError.resize(nPoint,nVar) = su2double(0.0); - /*--- Size Grad_AuxVar for axiysmmetric ---*/ if (config->GetAxisymmetric()){ nAuxVar = 3; @@ -83,29 +74,7 @@ CNEMOEulerVariable::CNEMOEulerVariable(su2double val_pressure, AuxVar.resize(nPoint,nAuxVar) = su2double(0.0); } - /*--- Only for residual smoothing (multigrid) ---*/ - for (unsigned long iMesh = 0; iMesh <= config->GetnMGLevels(); iMesh++) { - if (config->GetMG_CorrecSmooth(iMesh) > 0) { - Residual_Sum.resize(nPoint,nVar); - Residual_Old.resize(nPoint,nVar); - break; - } - } - - /*--- Allocate undivided laplacian (centered) and limiter (upwind)---*/ - if (config->GetKind_ConvNumScheme_Flow() == SPACE_CENTERED) - Undivided_Laplacian.resize(nPoint,nVar); - - /*--- Always allocate the slope limiter, - and the auxiliar variables (check the logic - JST with 2nd order Turb model - ) ---*/ - Limiter.resize(nPoint,nVar) = su2double(0.0); - Limiter_Primitive.resize(nPoint,nPrimVarGrad) = su2double(0.0); - - Solution_Max.resize(nPoint,nPrimVarGrad) = su2double(0.0); - Solution_Min.resize(nPoint,nPrimVarGrad) = su2double(0.0); - /*--- Primitive and secondary variables ---*/ - Primitive.resize(nPoint,nPrimVar) = su2double(0.0); Primitive_Aux.resize(nPoint,nPrimVar) = su2double(0.0); Secondary.resize(nPoint,nPrimVar) = su2double(0.0); @@ -116,32 +85,6 @@ CNEMOEulerVariable::CNEMOEulerVariable(su2double val_pressure, eves.resize(nPoint, nSpecies) = su2double(0.0); Gamma.resize(nPoint) = su2double(0.0); - /*--- Compressible flow, gradients primitive variables ---*/ - Gradient_Primitive.resize(nPoint,nPrimVarGrad,nDim,0.0); - Gradient.resize(nPoint,nVar,nDim,0.0); - - if (config->GetReconstructionGradientRequired()) { - Gradient_Aux.resize(nPoint,nPrimVarGrad,nDim,0.0); - } - - if (config->GetKind_Gradient_Method() == WEIGHTED_LEAST_SQUARES) { - Rmatrix.resize(nPoint,nDim,nDim,0.0); - } - - Velocity2.resize(nPoint) = su2double(0.0); - Max_Lambda_Inv.resize(nPoint) = su2double(0.0); - Delta_Time.resize(nPoint) = su2double(0.0); - Lambda.resize(nPoint) = su2double(0.0); - Sensor.resize(nPoint) = su2double(0.0); - - /* Non-physical point (first-order) initialization. */ - Non_Physical.resize(nPoint) = false; - Non_Physical_Counter.resize(nPoint) = 0; - - /* Under-relaxation parameter. */ - UnderRelaxation.resize(nPoint) = su2double(1.0); - LocalCFL.resize(nPoint) = su2double(0.0); - /*--- Set mixture state ---*/ fluidmodel->SetTDStatePTTv(val_pressure, val_massfrac, val_temperature, val_temperature_ve); @@ -150,7 +93,7 @@ CNEMOEulerVariable::CNEMOEulerVariable(su2double val_pressure, const su2double soundspeed = fluidmodel->ComputeSoundSpeed(); const su2double sqvel = GeometryToolbox::SquaredNorm(nDim, val_mach) * pow(soundspeed,2); const auto& energies = fluidmodel->ComputeMixtureEnergies(); - + /*--- Loop over all points --*/ for(unsigned long iPoint = 0; iPoint < nPoint; ++iPoint){ @@ -316,5 +259,3 @@ bool CNEMOEulerVariable::Cons2PrimVar(su2double *U, su2double *V, return nonPhys; } - -void CNEMOEulerVariable::SetSolution_New() { Solution_New = Solution; } diff --git a/SU2_CFD/src/variables/CNEMONSVariable.cpp b/SU2_CFD/src/variables/CNEMONSVariable.cpp index d99c71c8b86a..727f6cf1efc7 100644 --- a/SU2_CFD/src/variables/CNEMONSVariable.cpp +++ b/SU2_CFD/src/variables/CNEMONSVariable.cpp @@ -30,7 +30,7 @@ CNEMONSVariable::CNEMONSVariable(su2double val_pressure, const su2double *val_massfrac, - su2double *val_mach, + const su2double *val_mach, su2double val_temperature, su2double val_temperature_ve, unsigned long npoint, @@ -38,7 +38,7 @@ CNEMONSVariable::CNEMONSVariable(su2double val_pressure, unsigned long val_nvar, unsigned long val_nvarprim, unsigned long val_nvarprimgrad, - CConfig *config, + const CConfig *config, CNEMOGas *fluidmodel) : CNEMOEulerVariable(val_pressure, val_massfrac, val_mach, @@ -52,8 +52,6 @@ CNEMONSVariable::CNEMONSVariable(su2double val_pressure, config, fluidmodel) { - - Temperature_Ref = config->GetTemperature_Ref(); Viscosity_Ref = config->GetViscosity_Ref(); Viscosity_Inf = config->GetViscosity_FreeStreamND(); @@ -74,32 +72,7 @@ CNEMONSVariable::CNEMONSVariable(su2double val_pressure, Roe_Dissipation.resize(nPoint) = su2double(0.0); Vortex_Tilting.resize(nPoint) = su2double(0.0); Max_Lambda_Visc.resize(nPoint) = su2double(0.0); -} - -bool CNEMONSVariable::SetVorticity(void) { - - for (unsigned long iPoint=0; iPointGetModVel_FreeStream() / config->GetRefLength(); @@ -42,17 +42,18 @@ CNSVariable::CNSVariable(su2double density, const su2double *velocity, su2double Roe_Dissipation.resize(nPoint) = su2double(0.0); Vortex_Tilting.resize(nPoint) = su2double(0.0); Max_Lambda_Visc.resize(nPoint) = su2double(0.0); + } void CNSVariable::SetRoe_Dissipation_NTS(unsigned long iPoint, su2double val_delta, su2double val_const_DES){ - static const su2double cnu = pow(0.09, 1.5), - ch1 = 3.0, - ch2 = 1.0, - ch3 = 2.0, - sigma_max = 1.0; + const su2double cnu = pow(0.09, 1.5), + ch1 = 3.0, + ch2 = 1.0, + ch3 = 2.0, + sigma_max = 1.0; unsigned long iDim; su2double Omega, Omega_2 = 0.0, Baux, Gaux, Lturb, Kaux, Aaux; @@ -219,4 +220,3 @@ void CNSVariable::SetSecondaryVar(unsigned long iPoint, CFluidModel *FluidModel) SetdktdT_rho( iPoint, FluidModel->GetdktdT_rho() ); } - diff --git a/SU2_CFD/src/variables/CScalarVariable.cpp b/SU2_CFD/src/variables/CScalarVariable.cpp index 42f275ee325e..95ab8041600b 100644 --- a/SU2_CFD/src/variables/CScalarVariable.cpp +++ b/SU2_CFD/src/variables/CScalarVariable.cpp @@ -27,7 +27,7 @@ #include "../../include/variables/CScalarVariable.hpp" -CScalarVariable::CScalarVariable(unsigned long npoint, unsigned long ndim, unsigned long nvar, CConfig* config) +CScalarVariable::CScalarVariable(unsigned long npoint, unsigned long ndim, unsigned long nvar, const CConfig* config) : CVariable(npoint, ndim, nvar, config), Gradient_Reconstruction(config->GetReconstructionGradientRequired() ? Gradient_Aux : Gradient) { /*--- Gradient related fields ---*/ diff --git a/SU2_CFD/src/variables/CVariable.cpp b/SU2_CFD/src/variables/CVariable.cpp index 62d18a52002e..e1703a5f04c7 100644 --- a/SU2_CFD/src/variables/CVariable.cpp +++ b/SU2_CFD/src/variables/CVariable.cpp @@ -25,12 +25,10 @@ * License along with SU2. If not, see . */ - #include "../../include/variables/CVariable.hpp" #include "../../../Common/include/parallelization/omp_structure.hpp" - -CVariable::CVariable(unsigned long npoint, unsigned long nvar, CConfig *config) { +CVariable::CVariable(unsigned long npoint, unsigned long nvar, const CConfig *config) { /*--- Initialize the number of solution variables. This version of the constructor will be used primarily for converting the @@ -46,7 +44,8 @@ CVariable::CVariable(unsigned long npoint, unsigned long nvar, CConfig *config) } -CVariable::CVariable(unsigned long npoint, unsigned long ndim, unsigned long nvar, CConfig *config, bool adjoint) { +CVariable::CVariable(unsigned long npoint, unsigned long ndim, unsigned long nvar, + const CConfig *config, bool adjoint) { /*--- Initializate the number of dimension and number of variables ---*/ nPoint = npoint; diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index dbaab919878d..c1490f906fe1 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -71,7 +71,7 @@ def main(): invwedge.cfg_dir = "nonequilibrium/invwedge" invwedge.cfg_file = "invwedge.cfg" invwedge.test_iter = 10 - invwedge.test_vals = [-0.956173,-1.480936,-16.738781,-17.063703,-17.011887,2.371977,1.732488,5.399642,0.953492] + invwedge.test_vals = [-0.954417, -1.479835, -16.737319, -17.063704, -17.010424, 2.372843, 1.754751, 5.401691, 0.954721] invwedge.su2_exec = "mpirun -n 2 SU2_CFD" invwedge.timeout = 1600 invwedge.new_output = True @@ -83,7 +83,7 @@ def main(): visc_cone.cfg_dir = "nonequilibrium/axi_visccone" visc_cone.cfg_file = "axi_visccone.cfg" visc_cone.test_iter = 10 - visc_cone.test_vals = [-5.173078, -5.697841, -20.831296, -20.719164, -23.419769, -1.564084, -2.069040, 2.203924, -2.590729] + visc_cone.test_vals = [-5.175785, -5.700548, -20.705025, -20.592408, -22.851499, -1.540045, -2.104353, 2.197906, -2.576183] visc_cone.su2_exec = "mpirun -n 2 SU2_CFD" visc_cone.timeout = 1600 visc_cone.new_output = True diff --git a/TestCases/serial_regression.py b/TestCases/serial_regression.py index 7357ab104fd8..9778bc98620d 100644 --- a/TestCases/serial_regression.py +++ b/TestCases/serial_regression.py @@ -69,7 +69,7 @@ def main(): invwedge.cfg_dir = "nonequilibrium/invwedge" invwedge.cfg_file = "invwedge.cfg" invwedge.test_iter = 10 - invwedge.test_vals = [-0.956173,-1.480936,-16.738781,-17.063703,-17.011887,2.371977,1.732488, 5.399642,0.953492] + invwedge.test_vals = [-0.954415, -1.479832, -16.737319, -17.063704, -17.010424, 2.372845, 1.754743, 5.401692, 0.954723] invwedge.su2_exec = "SU2_CFD" invwedge.timeout = 1600 invwedge.new_output = True @@ -81,7 +81,7 @@ def main(): visc_cone.cfg_dir = "nonequilibrium/axi_visccone" visc_cone.cfg_file = "axi_visccone.cfg" visc_cone.test_iter = 10 - visc_cone.test_vals = [-5.173078, -5.697841, -20.831296, -20.719164, -23.419769, -1.564084, -2.069040, 2.203924, -2.590729] + visc_cone.test_vals = [-5.175741, -5.700504, -20.706185, -20.612638, -22.855655, -1.539882, -2.104390, 2.197875, -2.576097] visc_cone.su2_exec = "SU2_CFD" visc_cone.timeout = 1600 visc_cone.new_output = True