diff --git a/.github/workflows/regression.yml b/.github/workflows/regression.yml index 98cf7ed4bd80..6a8811b77b22 100644 --- a/.github/workflows/regression.yml +++ b/.github/workflows/regression.yml @@ -128,7 +128,7 @@ jobs: uses: docker://ghcr.io/su2code/su2/test-su2:220614-1237 with: # -t -c - args: -b ${{github.ref}} -t develop -c develop -s ${{matrix.testscript}} + args: -b ${{github.ref}} -t develop -c feature_cat_wall -s ${{matrix.testscript}} - name: Cleanup uses: docker://ghcr.io/su2code/su2/test-su2:220614-1237 with: diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 7c30b5168f2d..9e47dc85e567 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -311,7 +311,6 @@ class CConfig { su2double *Isothermal_Temperature; /*!< \brief Specified isothermal wall temperatures (static). */ su2double *HeatTransfer_Coeff; /*!< \brief Specified heat transfer coefficients. */ su2double *HeatTransfer_WallTemp; /*!< \brief Specified temperatures at infinity alongside heat transfer coefficients. */ - su2double *Wall_Catalycity; /*!< \brief Specified wall species mass-fractions for catalytic boundaries. */ su2double *Heat_Flux; /*!< \brief Specified wall heat fluxes. */ su2double *Roughness_Height; /*!< \brief Equivalent sand grain roughness for the marker according to config file. */ su2double *Displ_Value; /*!< \brief Specified displacement for displacement boundaries. */ @@ -1169,20 +1168,24 @@ class CConfig { unsigned short maxBasisDim, /*!< \brief Maximum number of POD basis dimensions. */ rom_save_freq; /*!< \brief Frequency of unsteady time steps to save. */ - unsigned short nSpecies, /*!< \brief Number of transported species equations (for NEMO and species transport)*/ + unsigned short nSpecies; /*!< \brief Number of transported species equations (for NEMO and species transport)*/ /* other NEMO configure options*/ - iWall_Catalytic, - nWall_Catalytic; /*!< \brief No of catalytic walls */ - su2double *Gas_Composition, /*!< \brief Initial mass fractions of flow [dimensionless] */ + unsigned short nSpecies_Cat_Wall, /*!< \brief No. of species for a catalytic wall. */ + iWall_Catalytic, /*!< \brief Iterator over catalytic walls. */ + nWall_Catalytic; /*!< \brief No. of catalytic walls. */ + su2double *Gas_Composition, /*!< \brief Initial mass fractions of flow [dimensionless]. */ + *Supercatalytic_Wall_Composition, /*!< \brief Supercatalytic wall mass fractions [dimensionless]. */ pnorm_heat; /*!< \brief pnorm for heat-flux. */ bool frozen, /*!< \brief Flag for determining if mixture is frozen. */ ionization, /*!< \brief Flag for determining if free electron gas is in the mixture. */ vt_transfer_res_limit, /*!< \brief Flag for determining if residual limiting for source term VT-transfer is used. */ - monoatomic; /*!< \brief Flag for monoatomic mixture. */ + monoatomic, /*!< \brief Flag for monoatomic mixture. */ + Supercatalytic_Wall; /*!< \brief Flag for supercatalytic wall. */ string GasModel, /*!< \brief Gas Model. */ *Wall_Catalytic; /*!< \brief Pointer to catalytic walls. */ TRANSCOEFFMODEL Kind_TransCoeffModel; /*!< \brief Transport coefficient Model for NEMO solver. */ + su2double CatalyticEfficiency; /*!< \brief Wall catalytic efficiency. */ /*--- Additional species solver options ---*/ bool Species_Clipping; /*!< \brief Boolean that activates solution clipping for scalar transport. */ @@ -3029,13 +3032,13 @@ class CConfig { /*! * \brief Retrieves the number of periodic time instances for Harmonic Balance. - * \return: Number of periodic time instances for Harmonic Balance. + * \return Number of periodic time instances for Harmonic Balance. */ unsigned short GetnTimeInstances(void) const { return nTimeInstances; } /*! * \brief Retrieves the period of oscillations to be used with Harmonic Balance. - * \return: Period for Harmonic Balance. + * \return Period for Harmonic Balance. */ su2double GetHarmonicBalance_Period(void) const { return HarmonicBalance_Period; } @@ -3724,6 +3727,12 @@ class CConfig { */ string GetWall_Catalytic_TagBound(unsigned short val_marker) const { return Wall_Catalytic[val_marker]; } + /*! + * \brief Get wall catalytic efficiency. + * \return wall catalytic efficiency value. + */ + su2double GetCatalytic_Efficiency(void) const { return CatalyticEfficiency; } + /*! * \brief Fluid model that we are using. * \return Fluid model that we are using. @@ -5211,22 +5220,22 @@ class CConfig { TIME_MARCHING GetTime_Marching() const { return TimeMarching; } /*! - * \brief Provides the number of species present in the plasma - * \return: The number of species present in the plasma, read from input file + * \brief Provides the number of species present in the gas mixture. + * \return The number of species present in the gas mixture. */ unsigned short GetnSpecies() const { return nSpecies; } - /*! - * \brief Get the wall heat flux on a constant heat flux boundary. - * \return The heat flux. + /*! + * \brief Provides the gas mass fractions of the flow. + * \return Gas Mass fractions. */ - const su2double *GetWall_Catalycity(void) const { return Wall_Catalycity; } + const su2double *GetGas_Composition(void) const { return Gas_Composition; } /*! - * \brief Provides the gas mass fractions of the flow - * \return: Gas Mass fractions + * \brief Provides the gas mass fractions at the wall for supercat wall. + * \return Supercat wall gas mass fractions. */ - const su2double *GetGas_Composition(void) const { return Gas_Composition; } + const su2double *GetSupercatalytic_Wall_Composition(void) const { return Supercatalytic_Wall_Composition; } /*! * \brief Provides the restart information. @@ -5304,6 +5313,11 @@ class CConfig { */ bool GetMonoatomic(void) const { return monoatomic; } + /*! + * \brief Indicates whether supercatalytic wall is used. + */ + bool GetSupercatalytic_Wall(void) const { return Supercatalytic_Wall; } + /*! * \brief Information about computing and plotting the equivalent area distribution. * \return TRUE or FALSE depending if we are computing the equivalent area. @@ -6217,17 +6231,30 @@ class CConfig { bool GetViscous_Wall(unsigned short iMarker) const; /*! - * \brief Determines if problem is adjoint - * \return true if Adjoint + * \brief Determines whether a marker with index iMarker is a catalytic boundary. + * \param iMarker + * \return it marker with index iMarker is a catalytic boundary. + */ + bool GetCatalytic_Wall(unsigned short iMarker) const; + + /*! + * \brief Determines if problem is adjoint. + * \return true if Adjoint. */ bool GetContinuous_Adjoint(void) const { return ContinuousAdjoint; } /*! - * \brief Determines if problem is viscous - * \return true if Viscous + * \brief Determines if problem is viscous. + * \return true if Viscous. */ bool GetViscous(void) const { return Viscous; } + /*! + * \brief Determines if problem has catalytic walls. + * \return true if catalytic walls are present. + */ + bool GetCatalytic(void) const { return nWall_Catalytic > 0; } + /*! * \brief Provides the index of the solution in the container. * \param[in] val_eqsystem - Equation that is being solved. diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 30f4c7810f1d..cc9ea692c76a 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -1057,7 +1057,6 @@ void CConfig::SetPointersNull(void) { Kind_TimeNumScheme = EULER_IMPLICIT; - Gas_Composition = nullptr; } void CConfig::SetConfig_Options() { @@ -1183,7 +1182,13 @@ void CConfig::SetConfig_Options() { addBoolOption("VT_RESIDUAL_LIMITING", vt_transfer_res_limit, false); /* DESCRIPTION: List of catalytic walls */ addStringListOption("CATALYTIC_WALL", nWall_Catalytic, Wall_Catalytic); - + /* DESCRIPTION: Specfify super-catalytic wall */ + addBoolOption("SUPERCATALYTIC_WALL", Supercatalytic_Wall, false); + /* DESCRIPTION: Wall mass fractions for supercatalytic case */ + addDoubleListOption("SUPERCATALYTIC_WALL_COMPOSITION", nSpecies_Cat_Wall, Supercatalytic_Wall_Composition); + /* DESCRIPTION: Specfify catalytic efficiency of wall if using gamma model */ + addDoubleOption("CATALYTIC_EFFICIENCY", CatalyticEfficiency, 1.0); + /*!\brief MARKER_MONITORING\n DESCRIPTION: Marker(s) of the surface where evaluate the non-dimensional coefficients \ingroup Config*/ /*--- Options related to VAN der WAALS MODEL and PENG ROBINSON ---*/ @@ -7796,6 +7801,17 @@ bool CConfig::GetViscous_Wall(unsigned short iMarker) const { Marker_All_KindBC[iMarker] == CHT_WALL_INTERFACE); } +bool CConfig::GetCatalytic_Wall(unsigned short iMarker) const { + + bool catalytic = false; + for (unsigned short iMarker_Catalytic = 0; iMarker_Catalytic < nWall_Catalytic; iMarker_Catalytic++){ + string Catalytic_Tag = Wall_Catalytic[iMarker_Catalytic]; + if (Catalytic_Tag == Marker_All_TagBound[iMarker]) { catalytic = true; break; } + } + + return catalytic; +} + bool CConfig::GetSolid_Wall(unsigned short iMarker) const { return GetViscous_Wall(iMarker) || diff --git a/SU2_CFD/include/fluid/CNEMOGas.hpp b/SU2_CFD/include/fluid/CNEMOGas.hpp index 222776e20a1b..3eb6828a6cb3 100644 --- a/SU2_CFD/include/fluid/CNEMOGas.hpp +++ b/SU2_CFD/include/fluid/CNEMOGas.hpp @@ -84,6 +84,8 @@ class CNEMOGas : public CFluidModel { Enthalpy_Formation, /*!< \brief Enthalpy of formation */ Ref_Temperature; /*!< \brief Reference temperature for thermodynamic relations */ + su2matrix CatRecombTable; /*!< \brief Table for catalytic wall recombination pairs. */ + public: /*! @@ -258,4 +260,10 @@ class CNEMOGas : public CFluidModel { * \brief Get species formation enthalpy. */ virtual const vector& GetSpeciesFormationEnthalpy() = 0; + + /*! + * \brief Get catalytic wall recombination indices and constants. + */ + inline const su2matrix& GetCatalyticRecombination() const {return CatRecombTable;} + }; diff --git a/SU2_CFD/include/solvers/CFVMFlowSolverBase.inl b/SU2_CFD/include/solvers/CFVMFlowSolverBase.inl index 1e84b6efae5f..e293e0286b79 100644 --- a/SU2_CFD/include/solvers/CFVMFlowSolverBase.inl +++ b/SU2_CFD/include/solvers/CFVMFlowSolverBase.inl @@ -2437,9 +2437,10 @@ void CFVMFlowSolverBase::Friction_Forces(const CGeometry* geometr unsigned long iVertex, iPoint, iPointNormal; unsigned short iMarker, iMarker_Monitoring, iDim, jDim; - su2double Viscosity = 0.0, Area, Density = 0.0, GradTemperature = 0.0, WallDistMod, FrictionVel, + su2double Viscosity = 0.0, Area, Density = 0.0, WallDistMod, FrictionVel, UnitNormal[3] = {0.0}, TauElem[3] = {0.0}, Tau[3][3] = {{0.0}}, Cp, - thermal_conductivity, MaxNorm = 8.0, Grad_Vel[3][3] = {{0.0}}, Grad_Temp[3] = {0.0}, AxiFactor; + thermal_conductivity, MaxNorm = 8.0, Grad_Vel[3][3] = {{0.0}}, Grad_Temp[3] = {0.0}, + Grad_Temp_ve[3] = {0.0}, AxiFactor; const su2double *Coord = nullptr, *Coord_Normal = nullptr, *Normal = nullptr; const su2double minYPlus = config->GetwallModel_MinYPlus(); @@ -2520,6 +2521,7 @@ void CFVMFlowSolverBase::Friction_Forces(const CGeometry* geometr Grad_Vel[iDim][jDim] = nodes->GetGradient_Primitive(iPoint, prim_idx.Velocity() + iDim, jDim); } Grad_Temp[iDim] = nodes->GetGradient_Primitive(iPoint, prim_idx.Temperature(), iDim); + if (nemo) Grad_Temp_ve[iDim] = nodes->GetGradient_Primitive(iPoint, prim_idx.Temperature_ve(), iDim); } Viscosity = nodes->GetLaminarViscosity(iPoint); @@ -2591,35 +2593,53 @@ void CFVMFlowSolverBase::Friction_Forces(const CGeometry* geometr /*--- Compute total and maximum heat flux on the wall ---*/ - /// TODO: Move these ifs to specialized functions. + su2double dTdn = -GeometryToolbox::DotProduct(nDim, Grad_Temp, UnitNormal); + if (!nemo){ if (FlowRegime == ENUM_REGIME::COMPRESSIBLE) { - GradTemperature = -GeometryToolbox::DotProduct(nDim, Grad_Temp, UnitNormal); Cp = (Gamma / Gamma_Minus_One) * Gas_Constant; thermal_conductivity = Cp * Viscosity / Prandtl_Lam; } if (FlowRegime == ENUM_REGIME::INCOMPRESSIBLE) { - if (energy) - GradTemperature = -GeometryToolbox::DotProduct(nDim, Grad_Temp, UnitNormal); - + if (!energy) dTdn = 0.0; thermal_conductivity = nodes->GetThermalConductivity(iPoint); } - HeatFlux[iMarker][iVertex] = -thermal_conductivity * GradTemperature * RefHeatFlux; + HeatFlux[iMarker][iVertex] = -thermal_conductivity * dTdn * RefHeatFlux; } else { const auto& thermal_conductivity_tr = nodes->GetThermalConductivity(iPoint); const auto& thermal_conductivity_ve = nodes->GetThermalConductivity_ve(iPoint); - const auto& Grad_PrimVar = nodes->GetGradient_Primitive(iPoint); - - su2double dTn = GeometryToolbox::DotProduct(nDim, Grad_PrimVar[prim_idx.Temperature()], UnitNormal); - su2double dTven = GeometryToolbox::DotProduct(nDim, Grad_PrimVar[prim_idx.Temperature_ve()], UnitNormal); + + const su2double dTvedn = -GeometryToolbox::DotProduct(nDim, Grad_Temp_ve, UnitNormal); + + /*--- Surface energy balance: trans-rot heat flux, vib-el heat flux ---*/ + HeatFlux[iMarker][iVertex] = -(thermal_conductivity_tr*dTdn + thermal_conductivity_ve*dTvedn); + + /*--- Compute enthalpy transport to surface due to mass diffusion ---*/ + bool catalytic = config->GetCatalytic_Wall(iMarker); + if (catalytic){ + + const auto nSpecies = config->GetnSpecies(); + const auto& Grad_PrimVar = nodes->GetGradient_Primitive(iPoint); + const auto& PrimVar = nodes->GetPrimitive(iPoint); + const auto& Ds = nodes->GetDiffusionCoeff(iPoint); + const auto& hs = nodes->GetEnthalpys(iPoint); + const su2double rho = PrimVar[prim_idx.Density()]; + + su2double sumJhs = 0.0; + for (auto iSpecies = 0u; iSpecies < nSpecies; iSpecies++) { + for (auto iDim = 0u; iDim < nDim; iDim++) { + su2double dYdn = 1.0/rho*(Grad_PrimVar[iSpecies][iDim] - PrimVar[iSpecies]*Grad_PrimVar[prim_idx.Density()][iDim]/rho); + sumJhs += rho*Ds[iSpecies]*hs[iSpecies]*dYdn*UnitNormal[iDim]; + } + } + /*--- Surface energy balance: mass diffusion ---*/ + HeatFlux[iMarker][iVertex] += sumJhs; - /*--- Surface energy balance: trans-rot heat flux, vib-el heat flux, - enthalpy transport due to mass diffusion ---*/ - HeatFlux[iMarker][iVertex] = thermal_conductivity_tr*dTn + thermal_conductivity_ve*dTven; + } } /*--- Note that heat is computed at the diff --git a/SU2_CFD/src/fluid/CMutationTCLib.cpp b/SU2_CFD/src/fluid/CMutationTCLib.cpp index 090d6b71fb47..32d252a55ef5 100644 --- a/SU2_CFD/src/fluid/CMutationTCLib.cpp +++ b/SU2_CFD/src/fluid/CMutationTCLib.cpp @@ -38,6 +38,7 @@ CMutationTCLib::CMutationTCLib(const CConfig* config, unsigned short val_nDim): Cv_ks.resize(nEnergyEq*nSpecies,0.0); es.resize(nEnergyEq*nSpecies,0.0); omega_vec.resize(1,0.0); + CatRecombTable.resize(nSpecies,2) = 0; /*--- Set up inputs to define type of mixture in the Mutation++ library ---*/ @@ -69,6 +70,33 @@ CMutationTCLib::CMutationTCLib(const CConfig* config, unsigned short val_nDim): } else { nHeavy = nSpecies; nEl = 0; } + /*--- Set up catalytic recombination table. ---*/ + // Creation/Destruction (+1/-1), Index of monoatomic reactants + // Monoatomic species (N,O) recombine into diaatomic (N2, O2) + if (gas_model == "N2") { + CatRecombTable(0,0) = 1; CatRecombTable(0,1) = 1; + CatRecombTable(1,0) = -1; CatRecombTable(1,1) = 1; + + } else if (gas_model == "air_5"){ + CatRecombTable(0,0) = -1; CatRecombTable(0,1) = 0; + CatRecombTable(1,0) = -1; CatRecombTable(1,1) = 1; + CatRecombTable(2,0) = 0; CatRecombTable(2,1) = 4; + CatRecombTable(3,0) = 1; CatRecombTable(3,1) = 0; + CatRecombTable(4,0) = 1; CatRecombTable(4,1) = 1; + + } else if (gas_model == "air_6") { + CatRecombTable(0,0) = -1; CatRecombTable(0,1) = 0; + CatRecombTable(1,0) = -1; CatRecombTable(1,1) = 1; + CatRecombTable(2,0) = 0; CatRecombTable(2,1) = 4; + CatRecombTable(3,0) = 1; CatRecombTable(3,1) = 0; + CatRecombTable(4,0) = 1; CatRecombTable(4,1) = 1; + CatRecombTable(5,0) = 0; CatRecombTable(5,1) = 4; + + } else { + if (config->GetCatalytic()) + SU2_MPI::Error("Cataylic wall recombination not implemented for specified Mutation gas model.", CURRENT_FUNCTION); + } + } CMutationTCLib::~CMutationTCLib(){} diff --git a/SU2_CFD/src/fluid/CSU2TCLib.cpp b/SU2_CFD/src/fluid/CSU2TCLib.cpp index e16f91b407c2..e65131958607 100644 --- a/SU2_CFD/src/fluid/CSU2TCLib.cpp +++ b/SU2_CFD/src/fluid/CSU2TCLib.cpp @@ -47,6 +47,7 @@ CSU2TCLib::CSU2TCLib(const CConfig* config, unsigned short val_nDim, bool viscou Omega00.resize(nSpecies,nSpecies,4,0.0); Omega11.resize(nSpecies,nSpecies,4,0.0); RxnConstantTable.resize(6,5) = su2double(0.0); + CatRecombTable.resize(nSpecies,2) = 0; Blottner.resize(nSpecies,3) = su2double(0.0); taus.resize(nSpecies,0.0); eve_eq.resize(nSpecies,0.0); @@ -110,6 +111,12 @@ CSU2TCLib::CSU2TCLib(const CConfig* config, unsigned short val_nDim, bool viscou ElDegeneracy(0,5) = 5; ElDegeneracy(0,6) = 15; + /*--- Catayltic wall table---*/ + // Creation/Destruction (+1/-1), Index of monoatomic reactants. + // Argon not used. + CatRecombTable(0,0) = 0; CatRecombTable(0,1) = 0; + + /*--- Values used in the Sutherland's formula. ---*/ if (viscous) { //F.M. White, Viscous Fluid Flow, 3rd ed., McGraw-Hill, 2006. mu_ref[0] = 2.125E-5; @@ -260,6 +267,13 @@ CSU2TCLib::CSU2TCLib(const CConfig* config, unsigned short val_nDim, bool viscou Omega11(1,0,0) = -8.3493693E-03; Omega11(1,0,1) = 1.7808911E-01; Omega11(1,0,2) = -1.4466155E+00; Omega11(1,0,3) = 1.9324210E+03; Omega11(1,1,0) = -7.7439615E-03; Omega11(1,1,1) = 1.7129007E-01; Omega11(1,1,2) = -1.4809088E+00; Omega11(1,1,3) = 2.1284951E+03; + /*--- Catayltic wall table---*/ + // Creation/Destruction (+1/-1), Index of monoatomic reactants. + // Monoatomic species (N,O) recombine into diaatomic (N2, O2) + CatRecombTable(0,0) = 1; CatRecombTable(0,1) = 1; + CatRecombTable(1,0) = -1; CatRecombTable(1,1) = 1; + + /*--- Values used in the Sutherland's formula. ---*/ if (viscous) { //F.M. White, Viscous Fluid Flow, 3rd ed., McGraw-Hill, 2006. k_ref[0] = 0.0242; @@ -614,6 +628,15 @@ CSU2TCLib::CSU2TCLib(const CConfig* config, unsigned short val_nDim, bool viscou Omega11(4,3,0) = -5.0478143E-03; Omega11(4,3,1) = 1.0236186E-01; Omega11(4,3,2) = -9.0058935E-01; Omega11(4,3,3) = 4.4472565E+02; Omega11(4,4,0) = -4.2451096E-03; Omega11(4,4,1) = 9.6820337E-02; Omega11(4,4,2) = -9.9770795E-01; Omega11(4,4,3) = 8.3320644E+02; + // Creation/Destruction (+1/-1), Index of monoatomic reactants + // Monoatomic species (N,O) recombine into diaatomic (N2, O2) + CatRecombTable(0,0) = 1; CatRecombTable(0,1) = 3; + CatRecombTable(1,0) = 1; CatRecombTable(1,1) = 4; + CatRecombTable(2,0) = 0; CatRecombTable(2,1) = 0; + CatRecombTable(3,0) = -1; CatRecombTable(3,1) = 3; + CatRecombTable(4,0) = -1; CatRecombTable(4,1) = 4; + + /*--- Values used in the Sutherland's formula. ---*/ if (viscous) { //F.M. White, Viscous Fluid Flow, 3rd ed., McGraw-Hill, 2006. k_ref[0] = 0.0241; @@ -1052,6 +1075,17 @@ CSU2TCLib::CSU2TCLib(const CConfig* config, unsigned short val_nDim, bool viscou Omega11(4,3,0) = -5.0478143E-03; Omega11(4,3,1) = 1.0236186E-01; Omega11(4,3,2) = -9.0058935E-01; Omega11(4,3,3) = 4.4472565E+02; Omega11(4,4,0) = -4.2451096E-03; Omega11(4,4,1) = 9.6820337E-02; Omega11(4,4,2) = -9.9770795E-01; Omega11(4,4,3) = 8.3320644E+02; + // Creation/Destruction (+1/-1), Index of monoatomic reactants + // Monoatomic species (N,O) recombine into diaatomic (N2, O2) + CatRecombTable(0,0) = 0; CatRecombTable(0,1) = 1; + CatRecombTable(1,0) = 1; CatRecombTable(1,1) = 4; + CatRecombTable(2,0) = 1; CatRecombTable(2,1) = 5; + CatRecombTable(3,0) = 0; CatRecombTable(3,1) = 1; + CatRecombTable(4,0) = -1; CatRecombTable(4,1) = 4; + CatRecombTable(5,0) = -1; CatRecombTable(5,1) = 5; + CatRecombTable(6,0) = 0; CatRecombTable(6,1) = 1; + + /*--- Values for Sutherland's formula. ---*/ if (viscous) { //F.M. White, Viscous Fluid Flow, 3rd ed., McGraw-Hill, 2006. k_ref[0] = 0.0241; diff --git a/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp b/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp index dd1f25df7cac..a347597a4a3e 100644 --- a/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp @@ -1129,7 +1129,7 @@ void CNEMOEulerSolver::SetNondimensionalization(CConfig *config, unsigned short /*--- For inviscid flow, energy is calculated from the specified FreeStream quantities using the proper gas law. ---*/ - Energy_FreeStream = energies[0] + 0.5*sqvel; + Energy_FreeStream = energies[0] + 0.5*sqvel; } diff --git a/SU2_CFD/src/solvers/CNEMONSSolver.cpp b/SU2_CFD/src/solvers/CNEMONSSolver.cpp index 2e8db0a44592..48079e95bc08 100644 --- a/SU2_CFD/src/solvers/CNEMONSSolver.cpp +++ b/SU2_CFD/src/solvers/CNEMONSSolver.cpp @@ -331,32 +331,14 @@ void CNEMONSSolver::BC_HeatFlux_Wall(CGeometry *geometry, CConfig *config, unsigned short val_marker) { - bool catalytic = false; - unsigned short iMarker_Catalytic = 0; - - const auto Marker_Tag = config->GetMarker_All_TagBound(val_marker); - - while( iMarker_Catalytic < config->GetnWall_Catalytic()){ - - string Catalytic_Tag = config->GetWall_Catalytic_TagBound(iMarker_Catalytic); + bool catalytic = config->GetCatalytic_Wall(val_marker); - if (Marker_Tag == Catalytic_Tag){ + if (catalytic) { BC_HeatFluxCatalytic_Wall(geometry, solver_container, conv_numerics, + sour_numerics, config, val_marker); - catalytic = true; - BC_HeatFluxCatalytic_Wall(geometry, solver_container, conv_numerics, - sour_numerics, config, val_marker); - break; - - } else { - - iMarker_Catalytic++; - - } + } else { BC_HeatFluxNonCatalytic_Wall(geometry, solver_container, conv_numerics, + sour_numerics, config, val_marker); } - - if(!catalytic) BC_HeatFluxNonCatalytic_Wall(geometry, solver_container, conv_numerics, - sour_numerics, config, val_marker); - } void CNEMONSSolver::BC_HeatFluxCatalytic_Wall(CGeometry *geometry, @@ -532,28 +514,14 @@ void CNEMONSSolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_con CNumerics *conv_numerics, CNumerics *sour_numerics, CConfig *config, unsigned short val_marker) { - bool catalytic = false; - unsigned short iMarker_Catalytic = 0; - - const auto Marker_Tag = config->GetMarker_All_TagBound(val_marker); + bool catalytic = config->GetCatalytic_Wall(val_marker); - while( iMarker_Catalytic < config->GetnWall_Catalytic()){ + if (catalytic) { BC_IsothermalCatalytic_Wall(geometry, solver_container, conv_numerics, + sour_numerics, config, val_marker); - string Catalytic_Tag = config->GetWall_Catalytic_TagBound(iMarker_Catalytic); - - if (Marker_Tag == Catalytic_Tag){ - - catalytic = true; - BC_IsothermalCatalytic_Wall(geometry, solver_container, conv_numerics, - sour_numerics, config, val_marker); - break; - } else { - iMarker_Catalytic++; - } + } else { BC_IsothermalNonCatalytic_Wall(geometry, solver_container, conv_numerics, + sour_numerics, config, val_marker); } - - if(!catalytic) BC_IsothermalNonCatalytic_Wall(geometry, solver_container, conv_numerics, - sour_numerics, config, val_marker); } void CNEMONSSolver::BC_IsothermalNonCatalytic_Wall(CGeometry *geometry, @@ -688,35 +656,22 @@ void CNEMONSSolver::BC_IsothermalCatalytic_Wall(CGeometry *geometry, CConfig *config, unsigned short val_marker) { - SU2_MPI::Error("BC_ISOTHERMAL with catalytic wall: Not operational in NEMO.", CURRENT_FUNCTION); - /*--- Call standard isothermal BC to apply no-slip and energy b.c.'s ---*/ BC_IsothermalNonCatalytic_Wall(geometry, solver_container, conv_numerics, sour_numerics, config, val_marker); - ///////////// FINITE DIFFERENCE METHOD /////////////// /*--- Local variables ---*/ - bool implicit; - unsigned short iDim, iSpecies, jSpecies, iVar, jVar, kVar; - unsigned long iVertex, iPoint, jPoint; - su2double rho, *eves, *dTdU, *dTvedU, *Cvve, *Normal, Area, Ru, RuSI, - dij, *Di, *Vi, *Vj, *Yj, *dYdn, SdYdn, **GradY, **dVdU; - const su2double *Yst; - vector hs, Cvtrs; + unsigned short iSpecies, jSpecies, iVar, jVar, kVar; + su2double **GradY, **dVdU; /*--- Assign booleans ---*/ - implicit = (config->GetKind_TimeIntScheme() == EULER_IMPLICIT); - - /*--- Set "Proportional control" coefficient ---*/ - //su2double pcontrol = 0.6; - - /*--- Get species mass fractions at the wall ---*/ - Yst = config->GetWall_Catalycity(); + const bool implicit = (config->GetKind_TimeIntScheme() == EULER_IMPLICIT); + const bool Supercatalytic_Wall = config->GetSupercatalytic_Wall(); /*--- Get universal information ---*/ - RuSI = UNIVERSAL_GAS_CONSTANT; - Ru = 1000.0*RuSI; - auto& Ms = FluidModel->GetSpeciesMolarMass(); + const su2double RuSI = UNIVERSAL_GAS_CONSTANT; + const su2double Ru = 1000.0*RuSI; + const auto& Ms = FluidModel->GetSpeciesMolarMass(); /*--- Get the locations of the primitive variables ---*/ const unsigned short RHOS_INDEX = nodes->GetRhosIndex(); @@ -725,8 +680,6 @@ void CNEMONSSolver::BC_IsothermalCatalytic_Wall(CGeometry *geometry, const unsigned short TVE_INDEX = nodes->GetTveIndex(); /*--- Allocate arrays ---*/ - Yj = new su2double[nSpecies]; - dYdn = new su2double[nSpecies]; GradY = new su2double*[nSpecies]; for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) GradY[iSpecies] = new su2double[nDim]; @@ -735,133 +688,158 @@ void CNEMONSSolver::BC_IsothermalCatalytic_Wall(CGeometry *geometry, dVdU[iVar] = new su2double[nVar]; /*--- Loop over all of the vertices on this boundary marker ---*/ - for(iVertex = 0; iVertex < geometry->nVertex[val_marker]; iVertex++) { - iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); + for (auto iVertex = 0u; iVertex < geometry->nVertex[val_marker]; iVertex++) { + + const auto iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); /*--- Check if the node belongs to the domain (i.e, not a halo node) ---*/ if (geometry->nodes->GetDomain(iPoint)) { /*--- Compute closest normal neighbor ---*/ - jPoint = geometry->vertex[val_marker][iVertex]->GetNormal_Neighbor(); + const auto jPoint = geometry->vertex[val_marker][iVertex]->GetNormal_Neighbor(); /*--- Compute distance between wall & normal neighbor ---*/ - dij = 0.0; - for (iDim = 0; iDim < nDim; iDim++) { - dij += (geometry->nodes->GetCoord(jPoint, iDim) - - geometry->nodes->GetCoord(iPoint, iDim)) - * (geometry->nodes->GetCoord(jPoint, iDim) - - geometry->nodes->GetCoord(iPoint, iDim)); - } - dij = sqrt(dij); - + const auto Coord_i = geometry->nodes->GetCoord(iPoint); + const auto Coord_j = geometry->nodes->GetCoord(jPoint); + const su2double dij = GeometryToolbox::Distance(nDim, Coord_i, Coord_j); /*--- Compute dual-grid area and boundary normal ---*/ - Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); - Area = GeometryToolbox::Norm(nDim, Normal); + const auto Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); + const su2double Area = GeometryToolbox::Norm(nDim, Normal); /*--- Initialize the viscous residual to zero ---*/ - for (iVar = 0; iVar < nVar; iVar++) - Res_Visc[iVar] = 0.0; + for (iVar = 0; iVar < nVar; iVar++) Res_Visc[iVar] = 0.0; /*--- Get primitive information ---*/ - Vi = nodes->GetPrimitive(iPoint); - Vj = nodes->GetPrimitive(jPoint); - Di = nodes->GetDiffusionCoeff(iPoint); - eves = nodes->GetEve(iPoint); - hs = FluidModel->ComputeSpeciesEnthalpy(Vi[T_INDEX], Vi[TVE_INDEX], eves); - for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) - Yj[iSpecies] = Vj[RHOS_INDEX+iSpecies]/Vj[RHO_INDEX]; - rho = Vi[RHO_INDEX]; - dTdU = nodes->GetdTdU(iPoint); - dTvedU = nodes->GetdTvedU(iPoint); - - /*--- Calculate normal derivative of mass fraction ---*/ - for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) - dYdn[iSpecies] = (Yst[iSpecies]-Yj[iSpecies])/dij; - - /*--- Calculate supplementary quantities ---*/ - SdYdn = 0.0; - for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) - SdYdn += rho*Di[iSpecies]*dYdn[iSpecies]; - - for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { - Res_Visc[iSpecies] = -(-rho*Di[iSpecies]*dYdn[iSpecies] - +Yst[iSpecies]*SdYdn )*Area; - Res_Visc[nSpecies+nDim] += (Res_Visc[iSpecies]*hs[iSpecies] )*Area; - Res_Visc[nSpecies+nDim+1] += (Res_Visc[iSpecies]*eves[iSpecies])*Area; - } - - /*--- Viscous contribution to the residual at the wall ---*/ - LinSysRes.SubtractBlock(iPoint, Res_Visc); + const auto& Vi = nodes->GetPrimitive(iPoint); + const auto& Vj = nodes->GetPrimitive(jPoint); + const auto& Di = nodes->GetDiffusionCoeff(iPoint); + const auto& eves = nodes->GetEve(iPoint); + const auto& hs = FluidModel->ComputeSpeciesEnthalpy(Vi[T_INDEX], Vi[TVE_INDEX], eves); + const su2double rho = Vi[RHO_INDEX]; + const auto& dTdU = nodes->GetdTdU(iPoint); + const auto& dTvedU = nodes->GetdTvedU(iPoint); - if (implicit) { + if (Supercatalytic_Wall) { - /*--- Initialize the transformation matrix ---*/ - for (iVar = 0; iVar < nVar; iVar++) - for (jVar = 0; jVar < nVar; jVar++) { - dVdU[iVar][jVar] = 0.0; - Jacobian_j[iVar][jVar] = 0.0; - Jacobian_i[iVar][jVar] = 0.0; - } + const auto& Yst = config->GetSupercatalytic_Wall_Composition(); - /*--- Populate transformation matrix ---*/ - // dYsdrk + /*--- Calculate supplementary quantities ---*/ + su2double dYdn = 0.0, SdYdn = 0.0; for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { - for (jSpecies = 0; jSpecies < nSpecies; jSpecies++) - dVdU[iSpecies][jSpecies] += -1.0/rho*Yst[iSpecies]; - dVdU[iSpecies][iSpecies] += 1.0/rho; - } - for (iVar = 0; iVar < nVar; iVar++) { - dVdU[nSpecies+nDim][iVar] = dTdU[iVar]; - dVdU[nSpecies+nDim+1][iVar] = dTvedU[iVar]; + dYdn = (Yst[iSpecies]-Vj[RHOS_INDEX+iSpecies]/Vj[RHO_INDEX])/dij; + SdYdn += rho*Di[iSpecies]*dYdn; } - /*--- Calculate supplementary quantities ---*/ - Cvtrs = FluidModel->GetSpeciesCvTraRot(); - Cvve = nodes->GetCvve(iPoint); - - /*--- Take the primitive var. Jacobian & store in Jac. jj ---*/ - // Species mass fraction + /*--- Calculate species residual at wall ---*/ for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { - for (jSpecies = 0; jSpecies < nSpecies; jSpecies++) - Jacobian_j[iSpecies][jSpecies] += -Yst[iSpecies]*rho*Di[jSpecies]/dij; - Jacobian_j[iSpecies][iSpecies] += rho*Di[iSpecies]/dij - SdYdn; + dYdn = (Yst[iSpecies]-Vj[RHOS_INDEX+iSpecies]/Vj[RHO_INDEX])/dij; + Res_Visc[iSpecies] = -(-rho*Di[iSpecies]*dYdn+Yst[iSpecies]*SdYdn)*Area; } - // Temperature - for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { - for (jSpecies = 0; jSpecies < nSpecies; jSpecies++) { - Jacobian_j[nSpecies+nDim][iSpecies] += Jacobian_j[jSpecies][iSpecies]*hs[iSpecies]; + if (implicit) { + /*--- Initialize the transformation matrix ---*/ + for (iVar = 0; iVar < nVar; iVar++) + for (jVar = 0; jVar < nVar; jVar++) { + dVdU[iVar][jVar] = 0.0; + Jacobian_j[iVar][jVar] = 0.0; + Jacobian_i[iVar][jVar] = 0.0; + } + + /*--- Populate transformation matrix ---*/ + // dYsdrk + for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { + for (jSpecies = 0; jSpecies < nSpecies; jSpecies++) + dVdU[iSpecies][jSpecies] += -1.0/rho*Yst[iSpecies]; + dVdU[iSpecies][iSpecies] += 1.0/rho; + } + for (iVar = 0; iVar < nVar; iVar++) { + dVdU[nSpecies+nDim][iVar] = dTdU[iVar]; + dVdU[nSpecies+nDim+1][iVar] = dTvedU[iVar]; + } + + /*--- Calculate supplementary quantities ---*/ + const auto& Cvtrs = FluidModel->GetSpeciesCvTraRot(); + const auto& Cvve = nodes->GetCvve(iPoint); + + /*--- Take the primitive var. Jacobian & store in Jac. jj ---*/ + // Species mass fraction + for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { + for (jSpecies = 0; jSpecies < nSpecies; jSpecies++) + Jacobian_j[iSpecies][jSpecies] += -Yst[iSpecies]*rho*Di[jSpecies]/dij; + Jacobian_j[iSpecies][iSpecies] += rho*Di[iSpecies]/dij - SdYdn; + } + + // Temperature + for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { + for (jSpecies = 0; jSpecies < nSpecies; jSpecies++) { + Jacobian_j[nSpecies+nDim][iSpecies] += Jacobian_j[jSpecies][iSpecies]*hs[iSpecies]; + } + Jacobian_j[nSpecies+nDim][nSpecies+nDim] += Res_Visc[iSpecies]/Area*(Ru/Ms[iSpecies] + + Cvtrs[iSpecies] ); + Jacobian_j[nSpecies+nDim][nSpecies+nDim+1] += Res_Visc[iSpecies]/Area*Cvve[iSpecies]; + } + + // Vib.-El. Temperature + for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { + for (jSpecies = 0; jSpecies < nSpecies; jSpecies++) + Jacobian_j[nSpecies+nDim+1][iSpecies] += Jacobian_j[jSpecies][iSpecies]*eves[iSpecies]; + Jacobian_j[nSpecies+nDim+1][nSpecies+nDim+1] += Res_Visc[iSpecies]/Area*Cvve[iSpecies]; } - Jacobian_j[nSpecies+nDim][nSpecies+nDim] += Res_Visc[iSpecies]/Area*(Ru/Ms[iSpecies] + - Cvtrs[iSpecies] ); - Jacobian_j[nSpecies+nDim][nSpecies+nDim+1] += Res_Visc[iSpecies]/Area*Cvve[iSpecies]; + + /*--- Multiply by the transformation matrix and store in Jac. ii ---*/ + for (iVar = 0; iVar < nVar; iVar++) + for (jVar = 0; jVar < nVar; jVar++) + for (kVar = 0; kVar < nVar; kVar++) + Jacobian_i[iVar][jVar] += Jacobian_j[iVar][kVar]*dVdU[kVar][jVar]*Area; + + /*--- Apply to the linear system ---*/ + Jacobian.SubtractBlock(iPoint, iPoint, Jacobian_i); } - // Vib.-El. Temperature - for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { - for (jSpecies = 0; jSpecies < nSpecies; jSpecies++) - Jacobian_j[nSpecies+nDim+1][iSpecies] += Jacobian_j[jSpecies][iSpecies]*eves[iSpecies]; - Jacobian_j[nSpecies+nDim+1][nSpecies+nDim+1] += Res_Visc[iSpecies]/Area*Cvve[iSpecies]; + } else { + + if (implicit) { + SU2_MPI::Error("Implicit not currently available for partially catalytic wall.",CURRENT_FUNCTION); } - /*--- Multiply by the transformation matrix and store in Jac. ii ---*/ - for (iVar = 0; iVar < nVar; iVar++) - for (jVar = 0; jVar < nVar; jVar++) - for (kVar = 0; kVar < nVar; kVar++) - Jacobian_i[iVar][jVar] += Jacobian_j[iVar][kVar]*dVdU[kVar][jVar]*Area; + /*--- Identify the boundary ---*/ + string Marker_Tag = config->GetMarker_All_TagBound(val_marker); - /*--- Apply to the linear system ---*/ - Jacobian.SubtractBlock2Diag(iPoint, Jacobian_i); + /*--- Get isothermal wall temp ----*/ + const su2double Tw = config->GetIsothermal_Temperature(Marker_Tag); + + /*--- Get wall catalytic efficiency ----*/ + const su2double gam = config->GetCatalytic_Efficiency(); + + /*--- Get cataltyic reaction map ---*/ + const auto& RxnTable = FluidModel->GetCatalyticRecombination(); + + /*--- Common catalytic flux factor ---*/ + const su2double factor = gam*rho*sqrt(RuSI*Tw/2/PI_NUMBER)*Area; + + /*--- Compute catalytic recombination flux ---*/ + // Ref: 10.2514/6.2022-1636 + // ws = gam_s*Ys*rho_wall*sqrt(Ru*Tw/(2*Pi*M_combine)*Area + for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { + int Index = SU2_TYPE::Int(RxnTable(iSpecies,1)); + Res_Visc[iSpecies] = RxnTable(iSpecies,0)*factor*Vi[Index]/Vi[RHO_INDEX]*sqrt(1/Ms[Index]); + } + } + + for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) { + Res_Visc[nSpecies+nDim] += (Res_Visc[iSpecies]*hs[iSpecies])*Area; + Res_Visc[nSpecies+nDim+1] += (Res_Visc[iSpecies]*eves[iSpecies])*Area; } + + /*--- Viscous contribution to the residual at the wall ---*/ + LinSysRes.SubtractBlock(iPoint, Res_Visc); } } for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) delete [] GradY[iSpecies]; delete [] GradY; - delete [] dYdn; - delete [] Yj; for (iVar = 0; iVar < nVar; iVar++) delete [] dVdU[iVar]; delete [] dVdU; diff --git a/TestCases/nonequilibrium/axi_visccone/axi_visccone.cfg b/TestCases/nonequilibrium/viscous/axi_visccone.cfg similarity index 97% rename from TestCases/nonequilibrium/axi_visccone/axi_visccone.cfg rename to TestCases/nonequilibrium/viscous/axi_visccone.cfg index 088452fcbd10..bc976a28601b 100644 --- a/TestCases/nonequilibrium/axi_visccone/axi_visccone.cfg +++ b/TestCases/nonequilibrium/viscous/axi_visccone.cfg @@ -1,7 +1,7 @@ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % % SU2 configuration file % -% Case description: Mach 5 viscous flow over a 10deg wedge % +% Case description: Mach 5 viscous flow over a 10deg axisymmetric cone % % Author: C. Garbacz % % Institution: Strathclyde University % % File Version 7.4.0 "Blackbird" % diff --git a/TestCases/nonequilibrium/viscous/partial_cat.cfg b/TestCases/nonequilibrium/viscous/partial_cat.cfg new file mode 100644 index 000000000000..d848fe962d28 --- /dev/null +++ b/TestCases/nonequilibrium/viscous/partial_cat.cfg @@ -0,0 +1,84 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Mach 5 viscous flow over a 10deg wedge with partially % +% catalytic walls - gamma model, effeciency = 0.2 % +% Author: J. Needels % +% Institution: Stanford University % +% File Version 7.4.0 "Blackbird" % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +SOLVER= NEMO_NAVIER_STOKES +GAS_MODEL= AIR-5 +GAS_COMPOSITION= (0.77, 0.23, 0.0, 0.0, 0.0) +MATH_PROBLEM= DIRECT +RESTART_SOL= NO +AXISYMMETRIC= NO + +% ----------- COMPRESSIBLE AND INCOMPRESSIBLE FREE-STREAM DEFINITION ----------% +% +MACH_NUMBER= 5 +AOA= 0.0 +SIDESLIP_ANGLE= 0.0 +INIT_OPTION= TD_CONDITIONS +FREESTREAM_PRESSURE= 10.0 +FREESTREAM_TEMPERATURE= 288.15 +FREESTREAM_TEMPERATURE_VE= 288.15 +REYNOLDS_NUMBER= 10000 +KIND_TURB_MODEL= NONE + +% ---- NONEQUILIBRIUM GAS, IDEAL GAS, POLYTROPIC, VAN DER WAALS AND PENG ROBINSON CONSTANTS -------% +% +FLUID_MODEL= SU2_NONEQ + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +% Euler wall boundary marker(s) (NONE = no marker) +MARKER_EULER= ( Euler) +MARKER_ISOTHERMAL= (Wall, 300) +MARKER_OUTLET= ( Exit, 5 ) +MARKER_FAR= ( Farfield, Inlet ) + +CATALYTIC_WALL= (wall) +SUPERCATALYTIC_WALL= NO +CATALYTIC_EFFICIENCY= 0.2 + +MARKER_PLOTTING= (NONE ) +MARKER_MONITORING= ( Wall ) + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +CFL_NUMBER= 0.9 +ITER= 11 +LINEAR_SOLVER= BCGSTAB +LINEAR_SOLVER_ERROR= 1E-6 +LINEAR_SOLVER_ITER= 5 + +% -----------------------------------------------------------------------% +% +CONV_NUM_METHOD_FLOW= AUSM +MUSCL_FLOW= NO +TIME_DISCRE_FLOW= EULER_EXPLICIT + +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +CONV_RESIDUAL_MINVAL= -15 +CONV_STARTITER= 10 + +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +MESH_FILENAME= viscwedge.su2 +MESH_FORMAT= SU2 +SOLUTION_FILENAME= restart_flow.dat +TABULAR_FORMAT= TECPLOT +CONV_FILENAME= convergence +RESTART_FILENAME= restart_flow.dat +VOLUME_FILENAME= soln_volume +SURFACE_FILENAME= soln_surface +OUTPUT_WRT_FREQ= 100 +OUTPUT_FILES= (RESTART_ASCII, PARAVIEW_ASCII) +SCREEN_OUTPUT= (INNER_ITER, RMS_DENSITY_0, RMS_DENSITY_1, RMS_DENSITY_2, RMS_DENSITY_3, RMS_DENSITY_4, RMS_ENERGY, RMS_ENERGY_VE, LIFT, DRAG, TOTAL_HEATFLUX) diff --git a/TestCases/nonequilibrium/viscous/super_cat.cfg b/TestCases/nonequilibrium/viscous/super_cat.cfg new file mode 100644 index 000000000000..a4bde59a0ea4 --- /dev/null +++ b/TestCases/nonequilibrium/viscous/super_cat.cfg @@ -0,0 +1,84 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Mach 5 viscous flow over a 10deg wedge with partially % +% catalytic walls - gamma model, effeciency = 0.2 % +% Author: J. Needels % +% Institution: Stanford University % +% File Version 7.4.0 "Blackbird" % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +SOLVER= NEMO_NAVIER_STOKES +GAS_MODEL= AIR-5 +GAS_COMPOSITION= (0.77, 0.23, 0.0, 0.0, 0.0) +MATH_PROBLEM= DIRECT +RESTART_SOL= NO +AXISYMMETRIC= NO + +% ----------- COMPRESSIBLE AND INCOMPRESSIBLE FREE-STREAM DEFINITION ----------% +% +MACH_NUMBER= 5 +AOA= 0.0 +SIDESLIP_ANGLE= 0.0 +INIT_OPTION= TD_CONDITIONS +FREESTREAM_PRESSURE= 10.0 +FREESTREAM_TEMPERATURE= 288.15 +FREESTREAM_TEMPERATURE_VE= 288.15 +REYNOLDS_NUMBER= 10000 +KIND_TURB_MODEL= NONE + +% ---- NONEQUILIBRIUM GAS, IDEAL GAS, POLYTROPIC, VAN DER WAALS AND PENG ROBINSON CONSTANTS -------% +% +FLUID_MODEL= SU2_NONEQ + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +% Euler wall boundary marker(s) (NONE = no marker) +MARKER_EULER= ( Euler) +MARKER_ISOTHERMAL= (Wall, 300) +MARKER_OUTLET= ( Exit, 5 ) +MARKER_FAR= ( Farfield, Inlet ) + +CATALYTIC_WALL= (wall) +SUPERCATALYTIC_WALL= YES +SUPERCATALYTIC_WALL_COMPOSITION= (0.77, 0.23, 0.0, 0.0, 0.0) + +MARKER_PLOTTING= (NONE ) +MARKER_MONITORING= ( Wall ) + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +CFL_NUMBER= 3.0 +ITER= 11 +LINEAR_SOLVER= BCGSTAB +LINEAR_SOLVER_ERROR= 1E-6 +LINEAR_SOLVER_ITER= 5 + +% -----------------------------------------------------------------------% +% +CONV_NUM_METHOD_FLOW= AUSM +MUSCL_FLOW= NO +TIME_DISCRE_FLOW= EULER_IMPLICIT + +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +CONV_RESIDUAL_MINVAL= -15 +CONV_STARTITER= 10 + +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +MESH_FILENAME= viscwedge.su2 +MESH_FORMAT= SU2 +SOLUTION_FILENAME= restart_flow.dat +TABULAR_FORMAT= TECPLOT +CONV_FILENAME= convergence +RESTART_FILENAME= restart_flow.dat +VOLUME_FILENAME= soln_volume +SURFACE_FILENAME= soln_surface +OUTPUT_WRT_FREQ= 100 +OUTPUT_FILES= (RESTART_ASCII, PARAVIEW_ASCII) +SCREEN_OUTPUT= (INNER_ITER, RMS_DENSITY_0, RMS_DENSITY_1, RMS_DENSITY_2, RMS_DENSITY_3, RMS_DENSITY_4, RMS_ENERGY, RMS_ENERGY_VE, LIFT, DRAG, TOTAL_HEATFLUX) \ No newline at end of file diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index 4cbd1bad3329..87b77b8e87db 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -82,7 +82,7 @@ def main(): # Viscous single cone - axisymmetric visc_cone = TestCase('visc_cone') - visc_cone.cfg_dir = "nonequilibrium/axi_visccone" + visc_cone.cfg_dir = "nonequilibrium/viscous" visc_cone.cfg_file = "axi_visccone.cfg" visc_cone.test_iter = 10 visc_cone.test_vals = [-5.222278, -5.746529, -20.569425, -20.633787, -20.547644, 1.255759, -3.208374, -0.016010, 0.093459, 32633.000000] @@ -99,6 +99,29 @@ def main(): #viscwedge_mpp.new_output = True #test_list.append(viscwedge_mpp) + # Viscous single wedge - super catalytic walls + super_cat = TestCase('super_cat') + super_cat.cfg_dir = "nonequilibrium/viscous" + super_cat.cfg_file = "super_cat.cfg" + super_cat.test_iter = 10 + super_cat.test_vals = [-5.232590, -5.757884, -20.727046, -20.748136, -20.564044, 1.246889, -3.205235, -0.028406, 0.250857, 3.2459e+04] + super_cat.su2_exec = "mpirun -n 2 SU2_CFD" + super_cat.timeout = 1600 + super_cat.new_output = True + super_cat.tol = 0.00001 + test_list.append(super_cat) + + # Viscous single wedge - partially catalytic walls + partial_cat = TestCase('partial_cat') + partial_cat.cfg_dir = "nonequilibrium/viscous" + partial_cat.cfg_file = "partial_cat.cfg" + partial_cat.test_iter = 10 + partial_cat.test_vals = [-5.210300, -5.735063, -20.880374, -20.825890, -23.475263, 1.806281, -2.813924, -0.078469, 0.496017, 2.9021e+04] + partial_cat.su2_exec = "mpirun -n 2 SU2_CFD" + partial_cat.timeout = 1600 + partial_cat.new_output = True + partial_cat.tol = 0.00001 + test_list.append(partial_cat) ########################## ### Compressible Euler ### diff --git a/TestCases/serial_regression.py b/TestCases/serial_regression.py index 3a2717ef3904..9685116e8535 100644 --- a/TestCases/serial_regression.py +++ b/TestCases/serial_regression.py @@ -70,7 +70,7 @@ def main(): # Viscous single cone - axisymmetric visc_cone = TestCase('visc_cone') - visc_cone.cfg_dir = "nonequilibrium/axi_visccone" + visc_cone.cfg_dir = "nonequilibrium/viscous" visc_cone.cfg_file = "axi_visccone.cfg" visc_cone.test_iter = 10 visc_cone.test_vals = [-5.215229, -5.739368, -20.545048, -20.618699, -20.502531, 1.262784, -3.205454, -0.015696, 0.093207, 32656.000000]