From dc66724638fae6bcb78953a66f93679dc77b9d1f Mon Sep 17 00:00:00 2001 From: bigfooted Date: Tue, 11 Nov 2025 22:55:28 +0100 Subject: [PATCH 01/31] custom species wall BC --- Common/include/CConfig.hpp | 22 +++++ Common/include/option_structure.hpp | 12 +++ Common/include/option_structure.inl | 78 ++++++++++++++++ Common/src/CConfig.cpp | 23 +++++ SU2_CFD/include/drivers/CDriverBase.hpp | 13 +++ .../include/solvers/CFVMFlowSolverBase.hpp | 8 ++ SU2_CFD/include/solvers/CSolver.hpp | 10 +++ SU2_CFD/include/solvers/CSpeciesSolver.hpp | 63 +++++++++++++ SU2_CFD/src/solvers/CSpeciesSolver.cpp | 88 ++++++++++++++++++- 9 files changed, 316 insertions(+), 1 deletion(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 794622bb1a1e..1d1d10a3fe6b 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -199,6 +199,7 @@ class CConfig { nMarker_Inlet, /*!< \brief Number of inlet flow markers. */ nMarker_Inlet_Species, /*!< \brief Number of inlet species markers. */ nSpecies_per_Inlet, /*!< \brief Number of species defined per inlet markers. */ + nMarker_Wall_Species, /*!< \brief Number of wall species markers. */ nMarker_Inlet_Turb, /*!< \brief Number of inlet turbulent markers. */ nTurb_Properties, /*!< \brief Number of turbulent properties per inlet markers. */ nMarker_Riemann, /*!< \brief Number of Riemann flow markers. */ @@ -255,6 +256,7 @@ class CConfig { *Marker_ActDiskBemOutlet_Axis, /*!< \brief Actuator disk BEM outlet markers passed to MARKER_ACTDISK_BEM_AXIS. */ *Marker_Inlet, /*!< \brief Inlet flow markers. */ *Marker_Inlet_Species, /*!< \brief Inlet species markers. */ + *Marker_Wall_Species, /*!< \brief Wall species markers. */ *Marker_Inlet_Turb, /*!< \brief Inlet turbulent markers. */ *Marker_Riemann, /*!< \brief Riemann markers. */ *Marker_Giles, /*!< \brief Giles markers. */ @@ -294,6 +296,8 @@ class CConfig { su2double **Inlet_Velocity; /*!< \brief Specified flow velocity vectors for supersonic inlet boundaries. */ su2double **Inlet_SpeciesVal; /*!< \brief Specified species vector for inlet boundaries. */ su2double **Inlet_TurbVal; /*!< \brief Specified turbulent intensity and viscosity ratio for inlet boundaries. */ + unsigned short *Kind_Wall_Species; /*!< \brief Species boundary condition type for wall boundaries (FLUX or VALUE). */ + su2double *Wall_SpeciesVal; /*!< \brief Specified species flux or value for wall boundaries. */ su2double *EngineInflow_Target; /*!< \brief Specified fan face targets for nacelle boundaries. */ su2double *Inflow_Mach; /*!< \brief Specified fan face mach for nacelle boundaries. */ su2double *Inflow_Pressure; /*!< \brief Specified fan face pressure for nacelle boundaries. */ @@ -1391,6 +1395,10 @@ class CConfig { void addGilesOption(const string name, unsigned short & nMarker_Giles, string * & Marker_Giles, unsigned short* & option_field, const map & enum_map, su2double* & var1, su2double* & var2, su2double** & FlowDir, su2double* & relaxfactor1, su2double* & relaxfactor2); + template + void addWallSpeciesOption(const string name, unsigned short & nMarker_Wall_Species, string * & Marker_Wall_Species, unsigned short* & option_field, const map & enum_map, + su2double* & value); + void addExhaustOption(const string& name, unsigned short & nMarker_Exhaust, string * & Marker_Exhaust, su2double* & Ttotal, su2double* & Ptotal); @@ -7131,6 +7139,20 @@ class CConfig { */ const su2double* GetInlet_SpeciesVal(const string& val_index) const; + /*! + * \brief Get the species value at a wall boundary + * \param[in] val_marker - Marker tag corresponding to the wall boundary. + * \return The wall species value (flux or Dirichlet value). + */ + su2double GetWall_SpeciesVal(const string& val_marker) const; + + /*! + * \brief Get the species boundary condition type at a wall boundary + * \param[in] val_marker - Marker tag corresponding to the wall boundary. + * \return The wall species type (WALL_SPECIES_FLUX or WALL_SPECIES_VALUE). + */ + unsigned short GetWall_SpeciesType(const string& val_marker) const; + /*! * \brief Get the turbulent properties values at an inlet boundary * \param[in] val_index - Index corresponding to the inlet boundary. diff --git a/Common/include/option_structure.hpp b/Common/include/option_structure.hpp index bf99257b1387..d8a560ab0fc4 100644 --- a/Common/include/option_structure.hpp +++ b/Common/include/option_structure.hpp @@ -1835,6 +1835,18 @@ static const MapType Giles_Map = { MakePair("MASS_FLOW_OUTLET", MASS_FLOW_OUTLET) }; +/*! + * \brief Types of wall species boundary conditions. + */ +enum WALL_SPECIES_TYPE { + WALL_SPECIES_FLUX = 1, /*!< \brief Neumann flux boundary condition for wall species. */ + WALL_SPECIES_VALUE = 2 /*!< \brief Dirichlet value boundary condition for wall species. */ +}; +static const MapType Wall_Map = { + MakePair("FLUX", WALL_SPECIES_FLUX) + MakePair("VALUE", WALL_SPECIES_VALUE) +}; + /*! * \brief Types of mixing process for averaging quantities at the boundaries. */ diff --git a/Common/include/option_structure.inl b/Common/include/option_structure.inl index 5f4a6e4fab51..5881620ab2d9 100644 --- a/Common/include/option_structure.inl +++ b/Common/include/option_structure.inl @@ -1302,6 +1302,84 @@ class COptionRiemann : public COptionBase { } }; +template +class COptionWallSpecies : public COptionBase { + protected: + map m; + string name; // identifier for the option + unsigned short& size; + string*& marker; + unsigned short*& field; // Reference to the field name + su2double*& value; + + public: + COptionWallSpecies(string option_field_name, unsigned short& nMarker_Wall_Species, string*& Marker_Wall_Species, + unsigned short*& option_field, const map m, su2double*& value) + : size(nMarker_Wall_Species), marker(Marker_Wall_Species), field(option_field), value(value) { + this->name = option_field_name; + this->m = m; + } + ~COptionWallSpecies() override{}; + + string SetValue(const vector& option_value) override { + COptionBase::SetValue(option_value); + unsigned short totalVals = option_value.size(); + if ((totalVals == 1) && (option_value[0].compare("NONE") == 0)) { + this->size = 0; + this->marker = nullptr; + this->field = nullptr; + this->value = nullptr; + return ""; + } + + if (totalVals % 3 != 0) { + string newstring; + newstring.append(this->name); + newstring.append(": must have a number of entries divisible by 3"); + this->size = 0; + this->marker = nullptr; + this->field = nullptr; + this->value = nullptr; + return newstring; + } + + unsigned short nVals = totalVals / 3; + this->size = nVals; + this->marker = new string[nVals]; + this->field = new unsigned short[nVals]; + this->value = new su2double[nVals]; + + for (unsigned long i = 0; i < nVals; i++) { + this->marker[i].assign(option_value[3 * i]); + // Check to see if the enum value is in the map + if (this->m.find(option_value[3 * i + 1]) == m.end()) { + string str; + str.append(this->name); + str.append(": invalid option value "); + str.append(option_value[3 * i + 1]); + str.append(". Check current SU2 options in config_template.cfg."); + return str; + } + Tenum val = this->m[option_value[3 * i + 1]]; + this->field[i] = val; + + istringstream ss_value(option_value[3 * i + 2]); + if (!(ss_value >> this->value[i])) { + return badValue("WallSpecies", this->name); + } + } + + return ""; + } + + void SetDefault() override { + this->marker = nullptr; + this->field = nullptr; + this->value = nullptr; + this->size = 0; // There is no default value for list + } +}; + template class COptionGiles : public COptionBase { map m; diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 0a9cca0f63f7..f5afada50191 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -513,6 +513,15 @@ void CConfig::addRiemannOption(const string name, unsigned short & nMarker_Riema option_map.insert(pair(name, val)); } +template +void CConfig::addWallSpeciesOption(const string name, unsigned short & nMarker_Wall_Species, string * & Marker_Wall_Species, unsigned short* & option_field, const map & enum_map, + su2double* & value) { + assert(option_map.find(name) == option_map.end()); + all_options.insert(pair(name, true)); + COptionBase* val = new COptionWallSpecies(name, nMarker_Wall_Species, Marker_Wall_Species, option_field, enum_map, value); + option_map.insert(pair(name, val)); +} + template void CConfig::addGilesOption(const string name, unsigned short & nMarker_Giles, string * & Marker_Giles, unsigned short* & option_field, const map & enum_map, su2double* & var1, su2double* & var2, su2double** & FlowDir, su2double* & relaxfactor1, su2double* & relaxfactor2) { @@ -9231,6 +9240,20 @@ const su2double* CConfig::GetInlet_SpeciesVal(const string& val_marker) const { return Inlet_SpeciesVal[iMarker_Inlet_Species]; } +su2double CConfig::GetWall_SpeciesVal(const string& val_marker) const { + unsigned short iMarker_Wall_Species; + for (iMarker_Wall_Species = 0; iMarker_Wall_Species < nMarker_Wall_Species; iMarker_Wall_Species++) + if (Marker_Wall_Species[iMarker_Wall_Species] == val_marker) break; + return Wall_SpeciesVal[iMarker_Wall_Species]; +} + +unsigned short CConfig::GetWall_SpeciesType(const string& val_marker) const { + unsigned short iMarker_Wall_Species; + for (iMarker_Wall_Species = 0; iMarker_Wall_Species < nMarker_Wall_Species; iMarker_Wall_Species++) + if (Marker_Wall_Species[iMarker_Wall_Species] == val_marker) break; + return Kind_Wall_Species[iMarker_Wall_Species]; +} + const su2double* CConfig::GetInlet_TurbVal(const string& val_marker) const { /*--- If Turbulent Inlet is not provided for the marker, return free stream values. ---*/ for (auto iMarker = 0u; iMarker < nMarker_Inlet_Turb; iMarker++) { diff --git a/SU2_CFD/include/drivers/CDriverBase.hpp b/SU2_CFD/include/drivers/CDriverBase.hpp index 97772527cf3f..5a55307b8819 100644 --- a/SU2_CFD/include/drivers/CDriverBase.hpp +++ b/SU2_CFD/include/drivers/CDriverBase.hpp @@ -566,6 +566,19 @@ class CDriverBase { main_geometry->SetCustomBoundaryHeatFlux(iMarker, iVertex, WallHeatFlux); } + /*! + * \brief Set the wall normal scalar values at a vertex on a specified marker (MARKER_PYTHON_CUSTOM). + * \note This can be the input of a scalar transport equation. + * \param[in] iMarker - Marker identifier. + * \param[in] iVertex - Vertex identifier. + * \param[in] WallScalar - Value of the normal heat flux. + */ + + inline void SetMarkerCustomScalar(unsigned short iMarker, unsigned long iVertex, std::vector WallScalar) { + auto* solver = solver_container[selected_zone][INST_0][MESH_0][SPECIES_SOL]; + solver->SetCustomBoundaryScalar(iMarker, iVertex, WallScalar); + } + /*! * \brief Selects zone to be used for python driver operations. * \param[in] iZone - Zone identifier. diff --git a/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp b/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp index e6f10ff2eaf7..a6f6daa95b2e 100644 --- a/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp +++ b/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp @@ -2200,6 +2200,14 @@ class CFVMFlowSolverBase : public CSolver { else Inlet_FlowDir[val_marker][val_vertex][val_dim] = val_flowdir; } + + /*! + * \brief Set the value of the customized normal scalar values/flux at a specified vertex on a specified marker. + * \param[in] val_marker - Marker value + * \param[in] val_vertex - Boundary vertex value + */ + inline void SetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, + std::vector val_customBoundaryScalar) { } /*! * \brief Update the multi-grid structure for the customized boundary conditions. diff --git a/SU2_CFD/include/solvers/CSolver.hpp b/SU2_CFD/include/solvers/CSolver.hpp index f0b3a4c493d1..6e864d0d7d1c 100644 --- a/SU2_CFD/include/solvers/CSolver.hpp +++ b/SU2_CFD/include/solvers/CSolver.hpp @@ -2855,6 +2855,16 @@ class CSolver { inline virtual su2double GetInletFlowDir(unsigned short val_marker, unsigned long val_vertex, unsigned short val_dim) const { return 0; } + + + /*! + * \brief Set the value of the customized normal scalar values/flux at a specified vertex on a specified marker. + * \param[in] val_marker - Marker value + * \param[in] val_vertex - Boundary vertex value + */ + inline virtual void SetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, + std::vector val_customBoundaryScalar) { } + /*! * \brief A virtual member diff --git a/SU2_CFD/include/solvers/CSpeciesSolver.hpp b/SU2_CFD/include/solvers/CSpeciesSolver.hpp index 052c3c21ee9c..67aa9053574f 100644 --- a/SU2_CFD/include/solvers/CSpeciesSolver.hpp +++ b/SU2_CFD/include/solvers/CSpeciesSolver.hpp @@ -40,6 +40,8 @@ class CSpeciesSolver : public CScalarSolver { protected: unsigned short Inlet_Position; /*!< \brief Column index for scalar variables in inlet files. */ vector Inlet_SpeciesVars; /*!< \brief Species variables at inlet profiles. */ + vector Wall_SpeciesVars; /*!< \brief Species variables at profiles. */ + vector >CustomBoundaryScalar; public: /*! @@ -146,6 +148,43 @@ class CSpeciesSolver : public CScalarSolver { void BC_Outlet(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, CNumerics* visc_numerics, CConfig* config, unsigned short val_marker) final; + /*! + * \brief Impose the isothermal wall boundary condition. + * \param[in] geometry - Geometrical definition of the problem. + * \param[in] solver_container - Container vector with all the solutions. + * \param[in] conv_numerics - Description of the numerical method. + * \param[in] visc_numerics - Description of the numerical method. + * \param[in] config - Definition of the particular problem. + * \param[in] val_marker - Surface marker where the boundary condition is applied. + */ + void BC_Isothermal_Wall(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, + CNumerics* visc_numerics, CConfig* config, unsigned short val_marker) override; + + /*! + * \brief Impose the heat flux wall boundary condition. + * \param[in] geometry - Geometrical definition of the problem. + * \param[in] solver_container - Container vector with all the solutions. + * \param[in] conv_numerics - Description of the numerical method. + * \param[in] visc_numerics - Description of the numerical method. + * \param[in] config - Definition of the particular problem. + * \param[in] val_marker - Surface marker where the boundary condition is applied. + */ + void BC_HeatFlux_Wall(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, + CNumerics* visc_numerics, CConfig* config, unsigned short val_marker) override; + + /*! + * \brief Generic wall boundary condition implementation. + * \param[in] geometry - Geometrical definition of the problem. + * \param[in] solver_container - Container vector with all the solutions. + * \param[in] conv_numerics - Description of the numerical method. + * \param[in] visc_numerics - Description of the numerical method. + * \param[in] config - Definition of the particular problem. + * \param[in] val_marker - Surface marker where the boundary condition is applied. + * \param[in] cht_mode - Flag for conjugate heat transfer mode (default: false). + */ + void BC_Wall_Generic(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, + CNumerics* visc_numerics, CConfig* config, unsigned short val_marker, bool cht_mode = false); + /*--- Note that BC_Sym_Plane, BC_HeatFlux_Wall, BC_Isothermal_Wall are all treated as zero-flux BC for the * mass-factions, which can be fulfilled by no additional residual contribution on these nodes. * If a specified mass-fractions flux (like BC_HeatFlux_Wall) or a constant mass-fraction on the boundary @@ -195,4 +234,28 @@ class CSpeciesSolver : public CScalarSolver { geometry, solver_container, conv_numerics, visc_numerics, config); } + /*! + * \brief Set custom boundary scalar values from Python. + * \param[in] val_marker - Boundary marker index + * \param[in] val_vertex - Boundary vertex index + * \param[in] val_customBoundaryScalar - Vector of scalar values + */ + inline void SetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, + std::vector val_customBoundaryScalar) { + for (auto iVar = 0u; iVar < nVar; iVar++) { + CustomBoundaryScalar[val_marker](val_vertex, iVar) = val_customBoundaryScalar[iVar]; + } + } + + /*! + * \brief Get custom boundary scalar value for a specific variable. + * \param[in] val_marker - Boundary marker index + * \param[in] val_vertex - Boundary vertex index + * \param[in] iVar - Variable index + * \return Custom boundary scalar value + */ + inline su2double GetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, unsigned short iVar) const { + return CustomBoundaryScalar[val_marker](val_vertex, iVar); + } + }; diff --git a/SU2_CFD/src/solvers/CSpeciesSolver.cpp b/SU2_CFD/src/solvers/CSpeciesSolver.cpp index 3f25af7d1183..bd46afd5b072 100644 --- a/SU2_CFD/src/solvers/CSpeciesSolver.cpp +++ b/SU2_CFD/src/solvers/CSpeciesSolver.cpp @@ -27,6 +27,7 @@ #include "../../include/solvers/CSpeciesSolver.hpp" +#include "../../../Common/include/option_structure.hpp" #include "../../../Common/include/parallelization/omp_structure.hpp" #include "../../../Common/include/toolboxes/geometry_toolbox.hpp" #include "../../include/solvers/CScalarSolver.inl" @@ -110,6 +111,8 @@ void CSpeciesSolver::Initialize(CGeometry* geometry, CConfig* config, unsigned s nDim = geometry->GetnDim(); + AllocVectorOfMatrices( nVertex, nVar,CustomBoundaryScalar); + if (iMesh == MESH_0 || config->GetMGCycle() == FULLMG_CYCLE) { /*--- Define some auxiliary vector related with the residual ---*/ @@ -174,6 +177,16 @@ void CSpeciesSolver::Initialize(CGeometry* geometry, CConfig* config, unsigned s } } } + + Wall_SpeciesVars.resize(nMarker); + for (unsigned long iMarker = 0; iMarker < nMarker; iMarker++) { + Wall_SpeciesVars[iMarker].resize(nVertex[iMarker], nVar); + for (unsigned long iVertex = 0; iVertex < nVertex[iMarker]; ++iVertex) { + for (unsigned short iVar = 0; iVar < nVar; iVar++) { + Wall_SpeciesVars[iMarker](iVertex, iVar) = Solution_Inf[iVar]; + } + } + } } @@ -404,6 +417,79 @@ void CSpeciesSolver::BC_Inlet(CGeometry* geometry, CSolver** solver_container, C END_SU2_OMP_FOR } + + +void CSpeciesSolver::BC_Isothermal_Wall(CGeometry* geometry, CSolver** solver_container, + CNumerics* conv_numerics, CNumerics* visc_numerics, CConfig* config, + unsigned short val_marker) { + BC_Wall_Generic(geometry, solver_container, conv_numerics, visc_numerics, config, val_marker); +} + +void CSpeciesSolver::BC_HeatFlux_Wall(CGeometry* geometry, CSolver** solver_container, + CNumerics* conv_numerics, CNumerics* visc_numerics, CConfig* config, + unsigned short val_marker) { + BC_Wall_Generic(geometry, solver_container, conv_numerics, visc_numerics, config, val_marker,0); +} + +void CSpeciesSolver::BC_Wall_Generic(CGeometry* geometry, CSolver** solver_container, + CNumerics* conv_numerics, CNumerics* visc_numerics, + CConfig* config, unsigned short val_marker, bool cht_mode) { + const bool implicit = config->GetKind_TimeIntScheme_Flow() == EULER_IMPLICIT; + const bool py_custom = config->GetMarker_All_PyCustom(val_marker); + + string Marker_Tag = config->GetMarker_All_TagBound(val_marker); + + // Get wall species boundary condition type and value for this marker + su2double WallSpeciesValue = config->GetWall_SpeciesVal(Marker_Tag); + unsigned short wallspeciestype = config->GetWall_SpeciesType(Marker_Tag); + + SU2_OMP_FOR_DYN(OMP_MIN_SIZE) + for (auto iVertex = 0u; iVertex < geometry->nVertex[val_marker]; iVertex++) { + + const auto iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); + + if (!geometry->nodes->GetDomain(iPoint)) continue; + + const auto Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); + + su2double Area = GeometryToolbox::Norm(nDim, Normal); + + for (auto iVar = 0u; iVar < nVar; iVar++) { + + su2double WallSpecies = WallSpeciesValue; + + /*--- Get the scalar values from the python wrapper. ---*/ + if (py_custom) { + WallSpecies = GetCustomBoundaryScalar(val_marker, iVertex, iVar); + } + + switch(wallspeciestype){ + case WALL_SPECIES_FLUX: + //Flux Boundary condition + LinSysRes(iPoint, iVar) -= WallSpecies * Area; + + break; + case WALL_SPECIES_VALUE: + //Dirichlet Strong Boundary Condition + nodes->SetSolution(iPoint, iVar, WallSpecies); + nodes->SetSolution_Old(iPoint, iVar, WallSpecies); + + LinSysRes(iPoint, iVar) = 0.0; + + if (implicit) { + unsigned long total_index = iPoint * nVar + iVar; + Jacobian.DeleteValsRowi(total_index); + } + break; + } + } + } + END_SU2_OMP_FOR + +} + + + void CSpeciesSolver::SetInletAtVertex(const su2double *val_inlet, unsigned short iMarker, unsigned long iVertex) { @@ -579,4 +665,4 @@ void CSpeciesSolver::SetInitialCondition(CGeometry **geometry, CSolver ***solver const bool restart = config->GetRestart() || config->GetRestart_Flow(); PushSolutionBackInTime(TimeIter, restart, solver_container, geometry, config); -} \ No newline at end of file +} From 1d811038273c8a879a6a2c1e0a9c75cff9220c90 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Tue, 11 Nov 2025 23:08:49 +0100 Subject: [PATCH 02/31] default zero flux species wall BC --- Common/src/CConfig.cpp | 24 ++++++++++++++++-------- config_template.cfg | 7 +++++++ 2 files changed, 23 insertions(+), 8 deletions(-) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index f5afada50191..5a8fef28e5c4 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -9241,17 +9241,25 @@ const su2double* CConfig::GetInlet_SpeciesVal(const string& val_marker) const { } su2double CConfig::GetWall_SpeciesVal(const string& val_marker) const { - unsigned short iMarker_Wall_Species; - for (iMarker_Wall_Species = 0; iMarker_Wall_Species < nMarker_Wall_Species; iMarker_Wall_Species++) - if (Marker_Wall_Species[iMarker_Wall_Species] == val_marker) break; - return Wall_SpeciesVal[iMarker_Wall_Species]; + /*--- Search for the marker in the wall species list ---*/ + for (unsigned short iMarker_Wall_Species = 0; iMarker_Wall_Species < nMarker_Wall_Species; iMarker_Wall_Species++) { + if (Marker_Wall_Species[iMarker_Wall_Species] == val_marker) { + return Wall_SpeciesVal[iMarker_Wall_Species]; + } + } + /*--- If marker not found (MARKER_WALL_SPECIES=NONE), return zero flux ---*/ + return 0.0; } unsigned short CConfig::GetWall_SpeciesType(const string& val_marker) const { - unsigned short iMarker_Wall_Species; - for (iMarker_Wall_Species = 0; iMarker_Wall_Species < nMarker_Wall_Species; iMarker_Wall_Species++) - if (Marker_Wall_Species[iMarker_Wall_Species] == val_marker) break; - return Kind_Wall_Species[iMarker_Wall_Species]; + /*--- Search for the marker in the wall species list ---*/ + for (unsigned short iMarker_Wall_Species = 0; iMarker_Wall_Species < nMarker_Wall_Species; iMarker_Wall_Species++) { + if (Marker_Wall_Species[iMarker_Wall_Species] == val_marker) { + return Kind_Wall_Species[iMarker_Wall_Species]; + } + } + /*--- If marker not found (MARKER_WALL_SPECIES=NONE), return FLUX type (zero flux BC) ---*/ + return WALL_SPECIES_FLUX; } const su2double* CConfig::GetInlet_TurbVal(const string& val_marker) const { diff --git a/config_template.cfg b/config_template.cfg index aa7a1a97fc97..a03a9c6a9d26 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -1051,6 +1051,13 @@ MARKER_HEATTRANSFER= ( NONE ) % Format: ( marker name, constant wall temperature (K), ... ) MARKER_ISOTHERMAL= ( NONE ) % +% Wall species boundary condition marker(s) for species transport (NONE = no marker) +% Specify either a Neumann flux boundary condition (FLUX) or a Dirichlet value boundary condition (VALUE) +% Format: ( marker name, BC_TYPE, value, ... ) +% where BC_TYPE is either FLUX (Neumann) or VALUE (Dirichlet) +% Example: MARKER_WALL_SPECIES= (wall1, FLUX, 0.5, wall2, VALUE, 0.3) +% MARKER_WALL_SPECIES= ( NONE ) +% % Far-field boundary marker(s) (NONE = no marker) MARKER_FAR= ( farfield ) % From 747c7deca254d07aedcac6a725a0932ba351283f Mon Sep 17 00:00:00 2001 From: bigfooted Date: Wed, 12 Nov 2025 08:33:29 +0100 Subject: [PATCH 03/31] add keyword and testcase --- Common/src/CConfig.cpp | 4 + .../species2_primitiveVenturi_flux_value.cfg | 131 ++++++++++++++++++ 2 files changed, 135 insertions(+) create mode 100644 TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_flux_value.cfg diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 5a8fef28e5c4..1702bd5b01da 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -1627,6 +1627,10 @@ void CConfig::SetConfig_Options() { /*!\brief MARKER_RIEMANN \n DESCRIPTION: Riemann boundary marker(s) with the following formats, a unit vector. * \n OPTIONS: See \link Riemann_Map \endlink. The variables indicated by the option and the flow direction unit vector must be specified. \ingroup Config*/ addRiemannOption("MARKER_RIEMANN", nMarker_Riemann, Marker_Riemann, Kind_Data_Riemann, Riemann_Map, Riemann_Var1, Riemann_Var2, Riemann_FlowDir); + /*!\brief MARKER_WALL_SPECIES \n DESCRIPTION: Wall species boundary marker(s) with the following format: + * (marker_name, BC_TYPE, value, ...) where BC_TYPE is either FLUX (Neumann) or VALUE (Dirichlet). + * \n OPTIONS: See \link Wall_Map \endlink. \ingroup Config*/ + addWallSpeciesOption("MARKER_WALL_SPECIES", nMarker_Wall_Species, Marker_Wall_Species, Kind_Wall_Species, Wall_Map, Wall_SpeciesVal); /*!\brief MARKER_GILES \n DESCRIPTION: Giles boundary marker(s) with the following formats, a unit vector. */ /* \n OPTIONS: See \link Giles_Map \endlink. The variables indicated by the option and the flow direction unit vector must be specified. \ingroup Config*/ addGilesOption("MARKER_GILES", nMarker_Giles, Marker_Giles, Kind_Data_Giles, Giles_Map, Giles_Var1, Giles_Var2, Giles_FlowDir, RelaxFactorAverage, RelaxFactorFourier); diff --git a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_flux_value.cfg b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_flux_value.cfg new file mode 100644 index 000000000000..f601cefc4668 --- /dev/null +++ b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_flux_value.cfg @@ -0,0 +1,131 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Species mixing with 2 species, i.e. 1 transport equations % +% Author: T. Kattmann % +% Institution: Bosch Thermotechniek B.V. % +% Date: 2021/10/14 % +% File Version 8.3.0 "Harrier" % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +SOLVER= INC_RANS +KIND_TURB_MODEL= SST +% +% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% +% +INC_DENSITY_MODEL= CONSTANT +INC_DENSITY_INIT= 1.1766 +% +INC_VELOCITY_INIT= ( 1.00, 0.0, 0.0 ) +% +INC_ENERGY_EQUATION= YES +INC_TEMPERATURE_INIT= 300.0 +% +INC_NONDIM= DIMENSIONAL +% +% -------------------- FLUID PROPERTIES ------------------------------------- % +% +FLUID_MODEL= CONSTANT_DENSITY +% +CONDUCTIVITY_MODEL= CONSTANT_CONDUCTIVITY +THERMAL_CONDUCTIVITY_CONSTANT= 0.0357 +% +PRANDTL_LAM= 0.72 +TURBULENT_CONDUCTIVITY_MODEL= NONE +PRANDTL_TURB= 0.90 +% +VISCOSITY_MODEL= CONSTANT_VISCOSITY +MU_CONSTANT= 1.716E-5 +% +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +MARKER_HEATFLUX= ( wall, 0.0 ) +MARKER_SYM= ( axis ) +% +INC_INLET_TYPE= VELOCITY_INLET VELOCITY_INLET +MARKER_INLET= ( gas_inlet, 300, 1.0, 1.0, 0.0, 0.0,\ + air_axial_inlet, 300, 1.0, 0.0, -1.0, 0.0 ) +MARKER_INLET_SPECIES= (gas_inlet, 0.1,\ + air_axial_inlet, 0.0 ) +MARKER_WALL_SPECIES=(axis,FLUX,0.1,wall,VALUE,0.5) +% +INC_OUTLET_TYPE= PRESSURE_OUTLET +MARKER_OUTLET= ( outlet, 0.0) +% +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +% +CFL_NUMBER= 60 +CFL_REDUCTION_SPECIES= 1.0 +CFL_REDUCTION_TURB= 1.0 +% +ITER= 10 +% +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +LINEAR_SOLVER= FGMRES +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1E-8 +LINEAR_SOLVER_ITER= 5 +% +% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% +% +CONV_NUM_METHOD_FLOW= FDS +MUSCL_FLOW= YES +SLOPE_LIMITER_FLOW = NONE +TIME_DISCRE_FLOW= EULER_IMPLICIT +% +% -------------------- SCALAR TRANSPORT ---------------------------------------% +% +KIND_SCALAR_MODEL= SPECIES_TRANSPORT +DIFFUSIVITY_MODEL= CONSTANT_DIFFUSIVITY +DIFFUSIVITY_CONSTANT= 0.001 +% +CONV_NUM_METHOD_SPECIES= BOUNDED_SCALAR +MUSCL_SPECIES= NO +SLOPE_LIMITER_SPECIES = NONE +% +TIME_DISCRE_SPECIES= EULER_IMPLICIT +% +SPECIES_INIT= 1.0 +SPECIES_CLIPPING= YES +SPECIES_CLIPPING_MIN= 0.0 +SPECIES_CLIPPING_MAX= 1.0 +% +% -------------------- TURBULENT TRANSPORT ---------------------------------------% +% +CONV_NUM_METHOD_TURB= BOUNDED_SCALAR +MUSCL_TURB= NO +% +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +CONV_FIELD= RMS_PRESSURE, RMS_VELOCITY-X, RMS_VELOCITY-Y, RMS_TKE, RMS_SPECIES +CONV_RESIDUAL_MINVAL= -18 +CONV_STARTITER= 10 +% +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +MESH_FILENAME= primitiveVenturi.su2 +SCREEN_OUTPUT= INNER_ITER WALL_TIME \ + RMS_PRESSURE RMS_VELOCITY-X RMS_VELOCITY-Y RMS_TKE RMS_DISSIPATION RMS_SPECIES_0 \ + LINSOL_ITER LINSOL_RESIDUAL \ + LINSOL_ITER_TURB LINSOL_RESIDUAL_TURB \ + LINSOL_ITER_SPECIES LINSOL_RESIDUAL_SPECIES SURFACE_SPECIES_VARIANCE +SCREEN_WRT_FREQ_INNER= 10 +% +HISTORY_OUTPUT= RMS_RES FLOW_COEFF LINSOL SPECIES_COEFF SPECIES_COEFF_SURF +MARKER_ANALYZE= outlet gas_inlet air_axial_inlet +% +OUTPUT_FILES= RESTART_ASCII, PARAVIEW_MULTIBLOCK +VOLUME_OUTPUT= RESIDUAL, PRIMITIVE +OUTPUT_WRT_FREQ= 1000 +% +RESTART_SOL= NO +SOLUTION_FILENAME= solution +RESTART_FILENAME= restart +% +WRT_PERFORMANCE= YES From f1ab1a66ca35ad380cd05b0b52382c57e26c81c0 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Wed, 12 Nov 2025 09:29:01 +0100 Subject: [PATCH 04/31] switch to passivedouble --- SU2_CFD/include/drivers/CDriverBase.hpp | 2 +- SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp | 6 +++--- SU2_CFD/include/solvers/CSolver.hpp | 8 ++++---- SU2_CFD/include/solvers/CSpeciesSolver.hpp | 2 +- 4 files changed, 9 insertions(+), 9 deletions(-) diff --git a/SU2_CFD/include/drivers/CDriverBase.hpp b/SU2_CFD/include/drivers/CDriverBase.hpp index 5a55307b8819..33b75fad0138 100644 --- a/SU2_CFD/include/drivers/CDriverBase.hpp +++ b/SU2_CFD/include/drivers/CDriverBase.hpp @@ -574,7 +574,7 @@ class CDriverBase { * \param[in] WallScalar - Value of the normal heat flux. */ - inline void SetMarkerCustomScalar(unsigned short iMarker, unsigned long iVertex, std::vector WallScalar) { + inline void SetMarkerCustomScalar(unsigned short iMarker, unsigned long iVertex, vector WallScalar) { auto* solver = solver_container[selected_zone][INST_0][MESH_0][SPECIES_SOL]; solver->SetCustomBoundaryScalar(iMarker, iVertex, WallScalar); } diff --git a/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp b/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp index a6f6daa95b2e..6ac60a8e2754 100644 --- a/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp +++ b/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp @@ -2200,14 +2200,14 @@ class CFVMFlowSolverBase : public CSolver { else Inlet_FlowDir[val_marker][val_vertex][val_dim] = val_flowdir; } - + /*! * \brief Set the value of the customized normal scalar values/flux at a specified vertex on a specified marker. * \param[in] val_marker - Marker value * \param[in] val_vertex - Boundary vertex value */ - inline void SetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, - std::vector val_customBoundaryScalar) { } + inline void SetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, + vector val_customBoundaryScalar) { } /*! * \brief Update the multi-grid structure for the customized boundary conditions. diff --git a/SU2_CFD/include/solvers/CSolver.hpp b/SU2_CFD/include/solvers/CSolver.hpp index 6e864d0d7d1c..8f3f9c4136aa 100644 --- a/SU2_CFD/include/solvers/CSolver.hpp +++ b/SU2_CFD/include/solvers/CSolver.hpp @@ -2855,15 +2855,15 @@ class CSolver { inline virtual su2double GetInletFlowDir(unsigned short val_marker, unsigned long val_vertex, unsigned short val_dim) const { return 0; } - - + + /*! * \brief Set the value of the customized normal scalar values/flux at a specified vertex on a specified marker. * \param[in] val_marker - Marker value * \param[in] val_vertex - Boundary vertex value */ - inline virtual void SetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, - std::vector val_customBoundaryScalar) { } + inline virtual void SetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, + vector val_customBoundaryScalar) { } /*! diff --git a/SU2_CFD/include/solvers/CSpeciesSolver.hpp b/SU2_CFD/include/solvers/CSpeciesSolver.hpp index 67aa9053574f..8b15f6dfa9ad 100644 --- a/SU2_CFD/include/solvers/CSpeciesSolver.hpp +++ b/SU2_CFD/include/solvers/CSpeciesSolver.hpp @@ -241,7 +241,7 @@ class CSpeciesSolver : public CScalarSolver { * \param[in] val_customBoundaryScalar - Vector of scalar values */ inline void SetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, - std::vector val_customBoundaryScalar) { + vector val_customBoundaryScalar) { for (auto iVar = 0u; iVar < nVar; iVar++) { CustomBoundaryScalar[val_marker](val_vertex, iVar) = val_customBoundaryScalar[iVar]; } From d779979bd0e30ce850d40be6a7deba96a782fdf7 Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 12 Nov 2025 09:35:53 +0100 Subject: [PATCH 05/31] Potential fix for code scanning alert no. 5749: Resource not released in destructor Co-authored-by: Copilot Autofix powered by AI <62310815+github-advanced-security[bot]@users.noreply.github.com> --- Common/include/option_structure.inl | 15 ++++++++++++++- 1 file changed, 14 insertions(+), 1 deletion(-) diff --git a/Common/include/option_structure.inl b/Common/include/option_structure.inl index 5881620ab2d9..cce18828861c 100644 --- a/Common/include/option_structure.inl +++ b/Common/include/option_structure.inl @@ -1319,7 +1319,20 @@ class COptionWallSpecies : public COptionBase { this->name = option_field_name; this->m = m; } - ~COptionWallSpecies() override{}; + ~COptionWallSpecies() override { + if (marker) { + delete[] marker; + marker = nullptr; + } + if (field) { + delete[] field; + field = nullptr; + } + if (value) { + delete[] value; + value = nullptr; + } + } string SetValue(const vector& option_value) override { COptionBase::SetValue(option_value); From 38d926d322407fdc10d7014328260fcaf73a0aa8 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Wed, 12 Nov 2025 11:23:42 +0100 Subject: [PATCH 06/31] typo in config_template.cfg --- config_template.cfg | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/config_template.cfg b/config_template.cfg index a03a9c6a9d26..f1af6eb4b7a8 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -1056,7 +1056,7 @@ MARKER_ISOTHERMAL= ( NONE ) % Format: ( marker name, BC_TYPE, value, ... ) % where BC_TYPE is either FLUX (Neumann) or VALUE (Dirichlet) % Example: MARKER_WALL_SPECIES= (wall1, FLUX, 0.5, wall2, VALUE, 0.3) -% MARKER_WALL_SPECIES= ( NONE ) +MARKER_WALL_SPECIES= ( NONE ) % % Far-field boundary marker(s) (NONE = no marker) MARKER_FAR= ( farfield ) From 88e29228ed8b2c08354e3ff2189a54822dd84d6d Mon Sep 17 00:00:00 2001 From: bigfooted Date: Fri, 14 Nov 2025 08:39:51 +0100 Subject: [PATCH 07/31] support per-species flux or value --- Common/include/CConfig.hpp | 20 ++- Common/include/option_structure.inl | 130 +++++++++++++---- Common/src/CConfig.cpp | 25 ++-- SU2_CFD/src/solvers/CSpeciesSolver.cpp | 8 +- .../species3_primitiveVenturi_flux_value.cfg | 138 ++++++++++++++++++ config_template.cfg | 7 +- 6 files changed, 273 insertions(+), 55 deletions(-) create mode 100644 TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 1d1d10a3fe6b..369ca2c8b8fe 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -200,6 +200,7 @@ class CConfig { nMarker_Inlet_Species, /*!< \brief Number of inlet species markers. */ nSpecies_per_Inlet, /*!< \brief Number of species defined per inlet markers. */ nMarker_Wall_Species, /*!< \brief Number of wall species markers. */ + nSpecies_per_Wall, /*!< \brief Number of species defined per wall markers. */ nMarker_Inlet_Turb, /*!< \brief Number of inlet turbulent markers. */ nTurb_Properties, /*!< \brief Number of turbulent properties per inlet markers. */ nMarker_Riemann, /*!< \brief Number of Riemann flow markers. */ @@ -296,8 +297,8 @@ class CConfig { su2double **Inlet_Velocity; /*!< \brief Specified flow velocity vectors for supersonic inlet boundaries. */ su2double **Inlet_SpeciesVal; /*!< \brief Specified species vector for inlet boundaries. */ su2double **Inlet_TurbVal; /*!< \brief Specified turbulent intensity and viscosity ratio for inlet boundaries. */ - unsigned short *Kind_Wall_Species; /*!< \brief Species boundary condition type for wall boundaries (FLUX or VALUE). */ - su2double *Wall_SpeciesVal; /*!< \brief Specified species flux or value for wall boundaries. */ + unsigned short **Kind_Wall_Species; /*!< \brief Species boundary condition type for wall boundaries (FLUX or VALUE) per species. */ + su2double **Wall_SpeciesVal; /*!< \brief Specified species flux or value for wall boundaries per species. */ su2double *EngineInflow_Target; /*!< \brief Specified fan face targets for nacelle boundaries. */ su2double *Inflow_Mach; /*!< \brief Specified fan face mach for nacelle boundaries. */ su2double *Inflow_Pressure; /*!< \brief Specified fan face pressure for nacelle boundaries. */ @@ -1396,8 +1397,9 @@ class CConfig { su2double* & var1, su2double* & var2, su2double** & FlowDir, su2double* & relaxfactor1, su2double* & relaxfactor2); template - void addWallSpeciesOption(const string name, unsigned short & nMarker_Wall_Species, string * & Marker_Wall_Species, unsigned short* & option_field, const map & enum_map, - su2double* & value); + void addWallSpeciesOption(const string name, unsigned short & nMarker_Wall_Species, string * & Marker_Wall_Species, + unsigned short** & option_field, const map & enum_map, + su2double** & value, unsigned short & nSpecies_per_Wall); void addExhaustOption(const string& name, unsigned short & nMarker_Exhaust, string * & Marker_Exhaust, su2double* & Ttotal, su2double* & Ptotal); @@ -7140,18 +7142,20 @@ class CConfig { const su2double* GetInlet_SpeciesVal(const string& val_index) const; /*! - * \brief Get the species value at a wall boundary + * \brief Get the species value at a wall boundary for a specific species * \param[in] val_marker - Marker tag corresponding to the wall boundary. + * \param[in] iSpecies - Species index. * \return The wall species value (flux or Dirichlet value). */ - su2double GetWall_SpeciesVal(const string& val_marker) const; + su2double GetWall_SpeciesVal(const string& val_marker, unsigned short iSpecies) const; /*! - * \brief Get the species boundary condition type at a wall boundary + * \brief Get the species boundary condition type at a wall boundary for a specific species * \param[in] val_marker - Marker tag corresponding to the wall boundary. + * \param[in] iSpecies - Species index. * \return The wall species type (WALL_SPECIES_FLUX or WALL_SPECIES_VALUE). */ - unsigned short GetWall_SpeciesType(const string& val_marker) const; + unsigned short GetWall_SpeciesType(const string& val_marker, unsigned short iSpecies) const; /*! * \brief Get the turbulent properties values at an inlet boundary diff --git a/Common/include/option_structure.inl b/Common/include/option_structure.inl index cce18828861c..e7d9c6b4853b 100644 --- a/Common/include/option_structure.inl +++ b/Common/include/option_structure.inl @@ -1309,13 +1309,16 @@ class COptionWallSpecies : public COptionBase { string name; // identifier for the option unsigned short& size; string*& marker; - unsigned short*& field; // Reference to the field name - su2double*& value; + unsigned short**& field; // Reference to the field name (now 2D: marker x species) + su2double**& value; // Now 2D: marker x species + unsigned short& nSpecies_per_Wall; public: COptionWallSpecies(string option_field_name, unsigned short& nMarker_Wall_Species, string*& Marker_Wall_Species, - unsigned short*& option_field, const map m, su2double*& value) - : size(nMarker_Wall_Species), marker(Marker_Wall_Species), field(option_field), value(value) { + unsigned short**& option_field, const map m, su2double**& value, + unsigned short& nSpecies_per_Wall) + : size(nMarker_Wall_Species), marker(Marker_Wall_Species), field(option_field), value(value), + nSpecies_per_Wall(nSpecies_per_Wall) { this->name = option_field_name; this->m = m; } @@ -1325,10 +1328,16 @@ class COptionWallSpecies : public COptionBase { marker = nullptr; } if (field) { + for (unsigned short i = 0; i < size; i++) { + delete[] field[i]; + } delete[] field; field = nullptr; } if (value) { + for (unsigned short i = 0; i < size; i++) { + delete[] value[i]; + } delete[] value; value = nullptr; } @@ -1342,43 +1351,103 @@ class COptionWallSpecies : public COptionBase { this->marker = nullptr; this->field = nullptr; this->value = nullptr; + this->nSpecies_per_Wall = 0; return ""; } - if (totalVals % 3 != 0) { + /*--- Determine the number of markers and species per marker. + * Format: marker1, TYPE1, value1, TYPE2, value2, ..., marker2, TYPE1, value1, ... + * Each marker name starts with a letter, each TYPE is an enum string (starts with letter), + * and each value is numeric. Pattern: marker, (TYPE, value) x N ---*/ + + vector marker_indices; // Indices where markers start + vector species_counts; // Number of species per marker + + // Find all marker positions (strings starting with a letter that are not TYPE keywords) + for (unsigned short i = 0; i < totalVals; i++) { + if (isalpha(option_value[i][0])) { + // Check if this could be a TYPE keyword (i.e., is it in the enum map?) + if (this->m.find(option_value[i]) != m.end()) { + continue; // This is a TYPE keyword, not a marker + } + // This is a marker name + marker_indices.push_back(i); + } + } + + if (marker_indices.empty()) { string newstring; newstring.append(this->name); - newstring.append(": must have a number of entries divisible by 3"); - this->size = 0; - this->marker = nullptr; - this->field = nullptr; - this->value = nullptr; + newstring.append(": no valid markers found"); return newstring; } - unsigned short nVals = totalVals / 3; - this->size = nVals; - this->marker = new string[nVals]; - this->field = new unsigned short[nVals]; - this->value = new su2double[nVals]; + // Calculate number of species for each marker + for (unsigned short i = 0; i < marker_indices.size(); i++) { + unsigned short start_idx = marker_indices[i] + 1; // Start after marker name + unsigned short end_idx = (i + 1 < marker_indices.size()) ? marker_indices[i + 1] : totalVals; + unsigned short entries = end_idx - start_idx; - for (unsigned long i = 0; i < nVals; i++) { - this->marker[i].assign(option_value[3 * i]); - // Check to see if the enum value is in the map - if (this->m.find(option_value[3 * i + 1]) == m.end()) { - string str; - str.append(this->name); - str.append(": invalid option value "); - str.append(option_value[3 * i + 1]); - str.append(". Check current SU2 options in config_template.cfg."); - return str; + // Each species needs 2 entries: TYPE and value + if (entries % 2 != 0) { + string newstring; + newstring.append(this->name); + newstring.append(": each marker must have pairs of (TYPE, value) entries"); + return newstring; } - Tenum val = this->m[option_value[3 * i + 1]]; - this->field[i] = val; - istringstream ss_value(option_value[3 * i + 2]); - if (!(ss_value >> this->value[i])) { - return badValue("WallSpecies", this->name); + species_counts.push_back(entries / 2); + } + + // Check that all markers have the same number of species + this->nSpecies_per_Wall = species_counts[0]; + for (auto count : species_counts) { + if (count != this->nSpecies_per_Wall) { + string newstring; + newstring.append(this->name); + newstring.append(": all markers must specify the same number of species"); + return newstring; + } + } + + // Allocate arrays + this->size = marker_indices.size(); + this->marker = new string[this->size]; + this->field = new unsigned short*[this->size]; + this->value = new su2double*[this->size]; + + for (unsigned short i = 0; i < this->size; i++) { + this->field[i] = new unsigned short[this->nSpecies_per_Wall]; + this->value[i] = new su2double[this->nSpecies_per_Wall]; + } + + // Parse the values + for (unsigned short iMarker = 0; iMarker < this->size; iMarker++) { + unsigned short marker_idx = marker_indices[iMarker]; + this->marker[iMarker].assign(option_value[marker_idx]); + + // Parse species data for this marker + for (unsigned short iSpecies = 0; iSpecies < this->nSpecies_per_Wall; iSpecies++) { + unsigned short type_idx = marker_idx + 1 + 2 * iSpecies; + unsigned short val_idx = type_idx + 1; + + // Check TYPE keyword + if (this->m.find(option_value[type_idx]) == m.end()) { + string str; + str.append(this->name); + str.append(": invalid option value "); + str.append(option_value[type_idx]); + str.append(". Check current SU2 options in config_template.cfg."); + return str; + } + + Tenum val = this->m[option_value[type_idx]]; + this->field[iMarker][iSpecies] = val; + + istringstream ss_value(option_value[val_idx]); + if (!(ss_value >> this->value[iMarker][iSpecies])) { + return badValue("WallSpecies", this->name); + } } } @@ -1390,6 +1459,7 @@ class COptionWallSpecies : public COptionBase { this->field = nullptr; this->value = nullptr; this->size = 0; // There is no default value for list + this->nSpecies_per_Wall = 0; } }; diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 1702bd5b01da..41ed3964c9bb 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -514,11 +514,12 @@ void CConfig::addRiemannOption(const string name, unsigned short & nMarker_Riema } template -void CConfig::addWallSpeciesOption(const string name, unsigned short & nMarker_Wall_Species, string * & Marker_Wall_Species, unsigned short* & option_field, const map & enum_map, - su2double* & value) { +void CConfig::addWallSpeciesOption(const string name, unsigned short & nMarker_Wall_Species, string * & Marker_Wall_Species, + unsigned short** & option_field, const map & enum_map, + su2double** & value, unsigned short & nSpecies_per_Wall) { assert(option_map.find(name) == option_map.end()); all_options.insert(pair(name, true)); - COptionBase* val = new COptionWallSpecies(name, nMarker_Wall_Species, Marker_Wall_Species, option_field, enum_map, value); + COptionBase* val = new COptionWallSpecies(name, nMarker_Wall_Species, Marker_Wall_Species, option_field, enum_map, value, nSpecies_per_Wall); option_map.insert(pair(name, val)); } @@ -1628,9 +1629,10 @@ void CConfig::SetConfig_Options() { * \n OPTIONS: See \link Riemann_Map \endlink. The variables indicated by the option and the flow direction unit vector must be specified. \ingroup Config*/ addRiemannOption("MARKER_RIEMANN", nMarker_Riemann, Marker_Riemann, Kind_Data_Riemann, Riemann_Map, Riemann_Var1, Riemann_Var2, Riemann_FlowDir); /*!\brief MARKER_WALL_SPECIES \n DESCRIPTION: Wall species boundary marker(s) with the following format: - * (marker_name, BC_TYPE, value, ...) where BC_TYPE is either FLUX (Neumann) or VALUE (Dirichlet). + * (marker_name, BC_TYPE, value, BC_TYPE, value, ...) where BC_TYPE is either FLUX (Neumann) or VALUE (Dirichlet). + * Each marker must specify the same number of species (N species per marker). * \n OPTIONS: See \link Wall_Map \endlink. \ingroup Config*/ - addWallSpeciesOption("MARKER_WALL_SPECIES", nMarker_Wall_Species, Marker_Wall_Species, Kind_Wall_Species, Wall_Map, Wall_SpeciesVal); + addWallSpeciesOption("MARKER_WALL_SPECIES", nMarker_Wall_Species, Marker_Wall_Species, Kind_Wall_Species, Wall_Map, Wall_SpeciesVal, nSpecies_per_Wall); /*!\brief MARKER_GILES \n DESCRIPTION: Giles boundary marker(s) with the following formats, a unit vector. */ /* \n OPTIONS: See \link Giles_Map \endlink. The variables indicated by the option and the flow direction unit vector must be specified. \ingroup Config*/ addGilesOption("MARKER_GILES", nMarker_Giles, Marker_Giles, Kind_Data_Giles, Giles_Map, Giles_Var1, Giles_Var2, Giles_FlowDir, RelaxFactorAverage, RelaxFactorFourier); @@ -5720,6 +5722,8 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i nSpecies_options.insert(nSpecies_options.end(), {nSpecies_Clipping_Min, nSpecies_Clipping_Max}); if (nMarker_Inlet_Species > 0) nSpecies_options.push_back(nSpecies_per_Inlet); + if (nMarker_Wall_Species > 0) + nSpecies_options.push_back(nSpecies_per_Wall); // Add more options for size check here. /*--- nSpecies_Init is the master, but it simply checks for consistency. ---*/ @@ -5732,6 +5736,7 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i /*--- Check whether some variables (or their sums) are in physical bounds. [0,1] for species related quantities. ---*/ /*--- Note, only for species transport, not for flamelet model ---*/ + /* if (Kind_Species_Model == SPECIES_MODEL::SPECIES_TRANSPORT) { su2double Species_Init_Sum = 0.0; for (unsigned short iSpecies = 0; iSpecies < nSpecies; iSpecies++) { @@ -5749,7 +5754,7 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i checkScalarBounds(Inlet_SpeciesVal_Sum, "MARKER_INLET_SPECIES sum", 0.0, 1.0); } } - +*/ } // species transport checks /*--- Define some variables for flamelet model. ---*/ @@ -9244,22 +9249,22 @@ const su2double* CConfig::GetInlet_SpeciesVal(const string& val_marker) const { return Inlet_SpeciesVal[iMarker_Inlet_Species]; } -su2double CConfig::GetWall_SpeciesVal(const string& val_marker) const { +su2double CConfig::GetWall_SpeciesVal(const string& val_marker, unsigned short iSpecies) const { /*--- Search for the marker in the wall species list ---*/ for (unsigned short iMarker_Wall_Species = 0; iMarker_Wall_Species < nMarker_Wall_Species; iMarker_Wall_Species++) { if (Marker_Wall_Species[iMarker_Wall_Species] == val_marker) { - return Wall_SpeciesVal[iMarker_Wall_Species]; + return Wall_SpeciesVal[iMarker_Wall_Species][iSpecies]; } } /*--- If marker not found (MARKER_WALL_SPECIES=NONE), return zero flux ---*/ return 0.0; } -unsigned short CConfig::GetWall_SpeciesType(const string& val_marker) const { +unsigned short CConfig::GetWall_SpeciesType(const string& val_marker, unsigned short iSpecies) const { /*--- Search for the marker in the wall species list ---*/ for (unsigned short iMarker_Wall_Species = 0; iMarker_Wall_Species < nMarker_Wall_Species; iMarker_Wall_Species++) { if (Marker_Wall_Species[iMarker_Wall_Species] == val_marker) { - return Kind_Wall_Species[iMarker_Wall_Species]; + return Kind_Wall_Species[iMarker_Wall_Species][iSpecies]; } } /*--- If marker not found (MARKER_WALL_SPECIES=NONE), return FLUX type (zero flux BC) ---*/ diff --git a/SU2_CFD/src/solvers/CSpeciesSolver.cpp b/SU2_CFD/src/solvers/CSpeciesSolver.cpp index bd46afd5b072..62154c735d4c 100644 --- a/SU2_CFD/src/solvers/CSpeciesSolver.cpp +++ b/SU2_CFD/src/solvers/CSpeciesSolver.cpp @@ -439,10 +439,6 @@ void CSpeciesSolver::BC_Wall_Generic(CGeometry* geometry, CSolver** solver_conta string Marker_Tag = config->GetMarker_All_TagBound(val_marker); - // Get wall species boundary condition type and value for this marker - su2double WallSpeciesValue = config->GetWall_SpeciesVal(Marker_Tag); - unsigned short wallspeciestype = config->GetWall_SpeciesType(Marker_Tag); - SU2_OMP_FOR_DYN(OMP_MIN_SIZE) for (auto iVertex = 0u; iVertex < geometry->nVertex[val_marker]; iVertex++) { @@ -456,6 +452,10 @@ void CSpeciesSolver::BC_Wall_Generic(CGeometry* geometry, CSolver** solver_conta for (auto iVar = 0u; iVar < nVar; iVar++) { + // Get wall species boundary condition type and value for this marker and species + su2double WallSpeciesValue = config->GetWall_SpeciesVal(Marker_Tag, iVar); + unsigned short wallspeciestype = config->GetWall_SpeciesType(Marker_Tag, iVar); + su2double WallSpecies = WallSpeciesValue; /*--- Get the scalar values from the python wrapper. ---*/ diff --git a/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg b/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg new file mode 100644 index 000000000000..e914bf4b1aa2 --- /dev/null +++ b/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg @@ -0,0 +1,138 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Species mixing with 3 species, i.e. 2 transport equations % +% using inlet file for 2 inlets % +% Author: T. Kattmann % +% Institution: Bosch Thermotechniek B.V. % +% Date: 2021/10/14 % +% File Version 8.3.0 "Harrier" % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +SOLVER= INC_RANS +KIND_TURB_MODEL= SST +% +% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% +% +INC_DENSITY_MODEL= CONSTANT +INC_DENSITY_INIT= 1.1766 +% +INC_VELOCITY_INIT= ( 1.00, 0.0, 0.0 ) +% +INC_ENERGY_EQUATION= YES +INC_TEMPERATURE_INIT= 300.0 +% +INC_NONDIM= INITIAL_VALUES +% +% -------------------- FLUID PROPERTIES ------------------------------------- % +% +FLUID_MODEL= CONSTANT_DENSITY +% +CONDUCTIVITY_MODEL= CONSTANT_CONDUCTIVITY +THERMAL_CONDUCTIVITY_CONSTANT= 0.0357 +% +PRANDTL_LAM= 0.72 +TURBULENT_CONDUCTIVITY_MODEL= NONE +PRANDTL_TURB= 0.90 +% +VISCOSITY_MODEL= CONSTANT_VISCOSITY +MU_CONSTANT= 1.716E-5 +% +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +MARKER_HEATFLUX= ( wall, 0.0 ) +MARKER_SYM= ( axis ) +% +SPECIFIED_INLET_PROFILE= NO +INLET_FILENAME= inlet_venturi.dat +INC_INLET_TYPE= VELOCITY_INLET VELOCITY_INLET +MARKER_INLET= ( gas_inlet, 300, 1.0, 1.0, 0.0, 0.0,\ + air_axial_inlet, 300, 1.0, 0.0, -1.0, 0.0 ) +MARKER_INLET_SPECIES= (gas_inlet, 0.5, 0.5,\ + air_axial_inlet, 0.6, 0.0 ) +MARKER_WALL_SPECIES=(axis,FLUX,0.1,FLUX,0.2,wall,VALUE,0.5,VALUE,0.6) +% +INC_OUTLET_TYPE= PRESSURE_OUTLET +MARKER_OUTLET= ( outlet, 0.0 ) +% +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +% +CFL_NUMBER= 2000 +CFL_REDUCTION_SPECIES= 1.0 +CFL_REDUCTION_TURB= 1.0 +% +% Run commented Iter for good results +%ITER= 1000 +ITER= 100 +% +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +LINEAR_SOLVER= FGMRES +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1E-8 +LINEAR_SOLVER_ITER= 5 +% +% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% +% +CONV_NUM_METHOD_FLOW= FDS +MUSCL_FLOW= YES +SLOPE_LIMITER_FLOW = NONE +TIME_DISCRE_FLOW= EULER_IMPLICIT +% +% -------------------- SCALAR TRANSPORT ---------------------------------------% +% +KIND_SCALAR_MODEL= SPECIES_TRANSPORT +DIFFUSIVITY_MODEL= CONSTANT_DIFFUSIVITY +DIFFUSIVITY_CONSTANT= 0.001 +% +CONV_NUM_METHOD_SPECIES= SCALAR_UPWIND +MUSCL_SPECIES= NO +SLOPE_LIMITER_SPECIES = NONE +% +TIME_DISCRE_SPECIES= EULER_IMPLICIT +% +SPECIES_INIT= 1.0, 0.0 +SPECIES_CLIPPING= YES +SPECIES_CLIPPING_MIN= 0.0, 0.0 +SPECIES_CLIPPING_MAX= 1.0, 1.0 +% +% -------------------- TURBULENT TRANSPORT ---------------------------------------% +% +CONV_NUM_METHOD_TURB= SCALAR_UPWIND +MUSCL_TURB= NO +% +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +CONV_FIELD= RMS_PRESSURE, RMS_VELOCITY-X, RMS_VELOCITY-Y, RMS_TKE, RMS_SPECIES +CONV_RESIDUAL_MINVAL= -18 +CONV_STARTITER= 10 +% +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +MESH_FILENAME= primitiveVenturi.su2 +SCREEN_OUTPUT= INNER_ITER WALL_TIME \ + RMS_PRESSURE RMS_VELOCITY-X RMS_VELOCITY-Y RMS_TKE RMS_DISSIPATION RMS_SPECIES_0 RMS_SPECIES_1 \ + LINSOL_ITER LINSOL_RESIDUAL \ + LINSOL_ITER_TURB LINSOL_RESIDUAL_TURB \ + LINSOL_ITER_SPECIES LINSOL_RESIDUAL_SPECIES +SCREEN_WRT_FREQ_INNER= 10 +% +HISTORY_OUTPUT= ITER RMS_RES LINSOL SPECIES_COEFF SPECIES_COEFF_SURF +CONV_FILENAME= history_inlet +MARKER_ANALYZE= gas_inlet, air_axial_inlet, outlet +MARKER_ANALYZE_AVERAGE= AREA +% +OUTPUT_FILES= RESTART_ASCII, PARAVIEW_MULTIBLOCK +VOLUME_OUTPUT= RESIDUAL, PRIMITIVE +OUTPUT_WRT_FREQ= 1000 +% +RESTART_SOL= NO +SOLUTION_FILENAME= solution +RESTART_FILENAME= restart +% +WRT_PERFORMANCE= YES diff --git a/config_template.cfg b/config_template.cfg index f1af6eb4b7a8..97f5d967970b 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -1053,9 +1053,10 @@ MARKER_ISOTHERMAL= ( NONE ) % % Wall species boundary condition marker(s) for species transport (NONE = no marker) % Specify either a Neumann flux boundary condition (FLUX) or a Dirichlet value boundary condition (VALUE) -% Format: ( marker name, BC_TYPE, value, ... ) -% where BC_TYPE is either FLUX (Neumann) or VALUE (Dirichlet) -% Example: MARKER_WALL_SPECIES= (wall1, FLUX, 0.5, wall2, VALUE, 0.3) +% for each species independently. All markers must specify the same number of species. +% Format: ( marker_name, BC_TYPE_1, value_1, BC_TYPE_2, value_2, ..., BC_TYPE_N, value_N, ... ) +% where BC_TYPE is either FLUX (Neumann) or VALUE (Dirichlet) for each species +% Example for 3 species: MARKER_WALL_SPECIES= (wall1, FLUX, 0.0, VALUE, 0.3, VALUE, 0.2, wall2, VALUE, 0.1, VALUE, 0.4, FLUX, 0.0) MARKER_WALL_SPECIES= ( NONE ) % % Far-field boundary marker(s) (NONE = no marker) From 19d0bd64b9dc1570d8d88b3709d52bd57bca65ad Mon Sep 17 00:00:00 2001 From: Nijso Date: Fri, 14 Nov 2025 16:05:33 +0100 Subject: [PATCH 08/31] Potential fix for code scanning alert no. 5752: Comparison of narrow type with wide type in loop condition Co-authored-by: Copilot Autofix powered by AI <62310815+github-advanced-security[bot]@users.noreply.github.com> --- Common/include/option_structure.inl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Common/include/option_structure.inl b/Common/include/option_structure.inl index e7d9c6b4853b..d3592564360e 100644 --- a/Common/include/option_structure.inl +++ b/Common/include/option_structure.inl @@ -1383,7 +1383,7 @@ class COptionWallSpecies : public COptionBase { } // Calculate number of species for each marker - for (unsigned short i = 0; i < marker_indices.size(); i++) { + for (size_t i = 0; i < marker_indices.size(); i++) { unsigned short start_idx = marker_indices[i] + 1; // Start after marker name unsigned short end_idx = (i + 1 < marker_indices.size()) ? marker_indices[i + 1] : totalVals; unsigned short entries = end_idx - start_idx; From 062e482ab290b3fbea3196a004af112fc6771a1d Mon Sep 17 00:00:00 2001 From: bigfooted Date: Fri, 14 Nov 2025 16:08:36 +0100 Subject: [PATCH 09/31] change out-of-bounds error into warning. --- Common/src/CConfig.cpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 41ed3964c9bb..d5e6d8b75f53 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -5711,7 +5711,7 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i /*--- Helper function that checks scalar variable bounds. ---*/ auto checkScalarBounds = [&](su2double scalar, const string& name, su2double lowerBound, su2double upperBound) { if (scalar < lowerBound || scalar > upperBound) - SU2_MPI::Error(string("Variable: ") + name + string(", is out of bounds."), CURRENT_FUNCTION); + cout << "Value: " << scalar << " for " << name << " is out of bounds [" << lowerBound << "," << upperBound << "]." << endl; }; /*--- Some options have to provide as many entries as there are additional species equations. ---*/ @@ -5736,7 +5736,7 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i /*--- Check whether some variables (or their sums) are in physical bounds. [0,1] for species related quantities. ---*/ /*--- Note, only for species transport, not for flamelet model ---*/ - /* + if (Kind_Species_Model == SPECIES_MODEL::SPECIES_TRANSPORT) { su2double Species_Init_Sum = 0.0; for (unsigned short iSpecies = 0; iSpecies < nSpecies; iSpecies++) { @@ -5754,7 +5754,7 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i checkScalarBounds(Inlet_SpeciesVal_Sum, "MARKER_INLET_SPECIES sum", 0.0, 1.0); } } -*/ + } // species transport checks /*--- Define some variables for flamelet model. ---*/ From 5f7c4c9fa2a723f522b60a1b801066cb7b5c0e71 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Sat, 15 Nov 2025 11:53:19 +0100 Subject: [PATCH 10/31] add testcases --- TestCases/py_wrapper/psi_heatloss/psi.cfg | 176 +++++++ TestCases/py_wrapper/psi_heatloss/run.py | 444 ++++++++++++++++++ .../zone_0/symmetry.vtu-934412289-428597.lock | Bin 0 -> 8 bytes 3 files changed, 620 insertions(+) create mode 100755 TestCases/py_wrapper/psi_heatloss/psi.cfg create mode 100755 TestCases/py_wrapper/psi_heatloss/run.py create mode 100644 TestCases/py_wrapper/psi_heatloss/vol_solution/zone_0/symmetry.vtu-934412289-428597.lock diff --git a/TestCases/py_wrapper/psi_heatloss/psi.cfg b/TestCases/py_wrapper/psi_heatloss/psi.cfg new file mode 100755 index 000000000000..6d1cb64275ed --- /dev/null +++ b/TestCases/py_wrapper/psi_heatloss/psi.cfg @@ -0,0 +1,176 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Species mixing with 3 species, i.e. 2 transport equations % +% Author: T. Kattmann % +% Institution: Bosch Thermotechniek B.V. % +% Date: 2021/10/14 % +% File Version 7.5.1 "Blackbird" % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +SOLVER= INC_RANS +KIND_TURB_MODEL= SST +SST_OPTIONS= V2003m + +RESTART_SOL= YES + +MGLEVEL= 0 +% + +% temperature and density +FREESTREAM_TEMPERATURE = 673 +FREESTREAM_DENSITY = 2.50 +% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% +% +%INC_DENSITY_MODEL= CONSTANT +INC_DENSITY_MODEL= VARIABLE +INC_DENSITY_INIT= 2.50 +% +INC_VELOCITY_INIT= (40.00, 0.0, 0.0 ) +% +INC_ENERGY_EQUATION= YES +INC_TEMPERATURE_INIT= 673.0 +% +INC_NONDIM= DIMENSIONAL +% +% -------------------- FLUID PROPERTIES ------------------------------------- % +% +%FLUID_MODEL= CONSTANT_DENSITY +FLUID_MODEL= INC_IDEAL_GAS +% +CONDUCTIVITY_MODEL= CONSTANT_CONDUCTIVITY +THERMAL_CONDUCTIVITY_CONSTANT= 0.0357 +% +PRANDTL_LAM= 0.72 +TURBULENT_CONDUCTIVITY_MODEL= NONE +PRANDTL_TURB= 0.90 +% +VISCOSITY_MODEL= SUTHERLAND +MU_CONSTANT= 1.716E-5 +MU_REF = 1.716e-5 +MU_T_REF= 273.15 +SUTHERLAND_CONSTANT = 110.4 + +SPECIFIC_HEAT_CP = 1150 +% +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +MARKER_HEATFLUX= ( wall_pipe,0.0, wall_interior, 0.0,wall_top, 0,wall_side,0 ) + +% note, case is axisymmetric +MARKER_SYM= ( symmetry ) +AXISYMMETRIC= YES +SPECIFIED_INLET_PROFILE= NO +INLET_FILENAME= inlet.dat +INC_INLET_TYPE= VELOCITY_INLET +INC_INLET_DAMPING= 0.01 +MARKER_INLET= ( inlet, 673, 40.0, 1.0, 0.0, 0.0) +MARKER_INLET_TURBULENT = (inlet, 0.10, 15) +MARKER_INLET_SPECIES= (inlet, 0.0, 304721) +INLET_MATCHING_TOLERANCE = 1e-3 +MARKER_PYTHON_CUSTOM= (wall_top) + +MARKER_WALL_SPECIES= wall_side, FLUX,0,VALUE,304721, wall_pipe, FLUX,0 , VALUE,304721 ,wall_interior , FLUX,0 , VALUE,304721 ,wall_top, FLUX,0,VALUE,-250000 + +MARKER_SPECIES_STRONG_BC=wall_top,wall_side,wall_pipe,wall_interior +% +INC_OUTLET_TYPE= PRESSURE_OUTLET +INC_OUTLET_DAMPING= 0.01 +MARKER_OUTLET= ( outlet, 0.0 ) +% +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +% +%CFL_NUMBER= 25.0 +CFL_NUMBER= 10.0 +CFL_REDUCTION_SPECIES= 1.0 +CFL_REDUCTION_TURB= 1.0 +CFL_ADAPT= NO +CFL_ADAPT_PARAM= ( 0.95, 1.01, 1.0, 250, 1.0e-4, 0) +% +% Run commented Iter for good results +ITER= 1 +% +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +LINEAR_SOLVER= FGMRES +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1E-10 +LINEAR_SOLVER_ITER= 10 + +% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% +% +CONV_NUM_METHOD_FLOW= FDS +MUSCL_FLOW= NO +SLOPE_LIMITER_FLOW = NONE +TIME_DISCRE_FLOW= EULER_IMPLICIT +% +% -------------------- SCALAR TRANSPORT ---------------------------------------% +% +KIND_SCALAR_MODEL= SPECIES_TRANSPORT +%DIFFUSIVITY_MODEL= CONSTANT_SCHMIDT +DIFFUSIVITY_MODEL= CONSTANT_DIFFUSIVITY +SCHMIDT_NUMBER_LAMINAR= 0.5 +DIFFUSIVITY_CONSTANT= 0.000148 + +% according to the paper +SCHMIDT_NUMBER_TURBULENT= 0.7 +% +CONV_NUM_METHOD_SPECIES= BOUNDED_SCALAR +MUSCL_SPECIES= NO +SLOPE_LIMITER_SPECIES = NONE +% +TIME_DISCRE_SPECIES= EULER_IMPLICIT +% +SPECIES_INIT= 0.0, 200000 + +SPECIES_CLIPPING= NO +SPECIES_CLIPPING_MIN= 0.0, -1.5e+6 +SPECIES_CLIPPING_MAX= 1.0, 304721 + +% +% -------------------- TURBULENT TRANSPORT ---------------------------------------% +% +CONV_NUM_METHOD_TURB= BOUNDED_SCALAR +MUSCL_TURB= NO + +% +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +CONV_FIELD= RMS_PRESSURE, RMS_VELOCITY-X, RMS_VELOCITY-Y, RMS_TKE, RMS_SPECIES +CONV_RESIDUAL_MINVAL= -18 +CONV_STARTITER= 10 +% +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +MESH_FILENAME= psi.su2 +%MESH_FILENAME= psi_medium.su2 +%MESH_FILENAME= psi_fine.su2 +% +SCREEN_OUTPUT= INNER_ITER WALL_TIME \ + RMS_PRESSURE RMS_VELOCITY-X RMS_VELOCITY-Y RMS_TKE RMS_DISSIPATION RMS_SPECIES_0 RMS_SPECIES_1 \ + LINSOL_ITER LINSOL_RESIDUAL \ + LINSOL_ITER_TURB LINSOL_RESIDUAL_TURB \ + LINSOL_ITER_SPECIES LINSOL_RESIDUAL_SPECIES \ + SURFACE_SPECIES_0 +SCREEN_WRT_FREQ_INNER= 1 +% +HISTORY_OUTPUT= ITER RMS_RES LINSOL SPECIES_COEFF SPECIES_COEFF_SURF +CONV_FILENAME= history +MARKER_ANALYZE= gas_inlet, air_axial_inlet, outlet +MARKER_ANALYZE_AVERAGE= AREA +% +OUTPUT_FILES= RESTART_ASCII, PARAVIEW_MULTIBLOCK +VOLUME_OUTPUT= RESIDUAL, PRIMITIVE, RANK, SPECIES_UDS_0 +OUTPUT_WRT_FREQ= 100 +% +READ_BINARY_RESTART= NO +RESTART_FILENAME= restart +SOLUTION_FILENAME= solution +% +WRT_PERFORMANCE= YES +PYTHON_CUSTOM_SOURCE= YES diff --git a/TestCases/py_wrapper/psi_heatloss/run.py b/TestCases/py_wrapper/psi_heatloss/run.py new file mode 100755 index 000000000000..92511a31a99a --- /dev/null +++ b/TestCases/py_wrapper/psi_heatloss/run.py @@ -0,0 +1,444 @@ +#!/usr/bin/env python + +## \file run.py +# \brief turbulent premixed dump combustor simulation (PSI flame) +# phi=0.5, methane-air, U=40 m/s +# \version 8.1.0 "Harrier" +# +# SU2 Project Website: https://su2code.github.io +# +# The SU2 Project is maintained by the SU2 Foundation +# (http://su2foundation.org) +# +# Copyright 2012-2024, 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 . + +import sys +import pysu2 +import numpy as np +from mpi4py import MPI + +# unburnt temperature of the propane-air mixture +# flame temperature of the methane-air mixture (phi=0.5, P=5) +Tf = 1803 +Twall=1200 +Tu = 673.0 +Pu = 5.0 +phi = 0.5 +# unburnt density at P=5 +rho_u = 2.50 +# unburnt thermal conductivity of methane-air at phi=0.5 (P=5) +k_u = 0.0523 +# unburnt heat capacity of methane-air at phi=0.5 (P=1) +cp_u = 1311.0 + +NASA_COEFFICIENTS = { + 'CH4': { + 'low_coeffs': [5.14911468, -0.013662201, 4.91454E-05, -4.84247E-08, 1.66603E-11, -10246.5983, -4.63848842], + 'high_coeffs': [1.65326226, 0.01002631, -3.31661E-06, 5.36483E-10, -3.14697E-14, -10009.5936, 9.90506283] + }, + 'O2': { + 'low_coeffs': [3.78245636, -0.002996734, 9.8473E-06, -9.6813E-09, 3.24373E-12, -1063.94356, 3.65767573], + 'high_coeffs': [3.66096083, 0.000656366, -1.41149E-07, 2.05798E-11, -1.29913E-15, -1215.97725, 3.41536184] + }, + 'N2': { + 'low_coeffs': [3.53100528, -0.000123661, -5.02999E-07, 2.43531E-09, -1.40881E-12, -1046.97628, 2.96747038], + 'high_coeffs': [2.95257637, 0.0013969, -4.92632E-07, 7.8601E-11, -4.60755E-15, -923.948688, 5.87188762] + }, + 'CO2': { + 'low_coeffs': [2.356813, 0.00898413, -7.12206E-06, 2.4573E-09, -1.42885E-13, -48371.971, 9.9009035], + 'high_coeffs': [4.6365111, 0.002741457, -9.95898E-07, 1.60387E-10, -9.16199E-15, -49024.904, -1.9348955] + }, + 'H2O': { + 'low_coeffs': [4.1986352, -0.002036402, 6.52034E-06, -5.48793E-09, 1.77197E-12, -30293.726, -0.84900901], + 'high_coeffs': [2.6770389, 0.002973182, -7.73769E-07, 9.44335E-11, -4.269E-15, -29885.894, 6.88255] + }, + 'H2': { + 'low_coeffs': [2.34433112, 0.007980521, -1.94782E-05, 2.01572E-08, -7.37612E-12, -917.935173, 0.683010238], + 'high_coeffs': [2.93286575, 0.000826608, -1.46402E-07, 1.541E-11, -6.88805E-16, -813.065581, -1.02432865] + } +} + +# Temperature range boundaries +T_LOW = 300.0 # K +T_HIGH = 1000.0 # K + +def get_nasa_coefficients(species, T): + """ + Return the appropriate coefficients for the given temperature + """ + if species not in NASA_COEFFICIENTS: + raise ValueError(f"Species {species} not found in database") + + if T <= T_HIGH: + return NASA_COEFFICIENTS[species]['low_coeffs'] + else: + return NASA_COEFFICIENTS[species]['high_coeffs'] + + +def calculate_cp_nasa(T, coeffs): + """ + Calculate Cp in kJ/kmol-K using NASA polynomial coefficients + Cp/R = a1 + a2*T + a3*T² + a4*T³ + a5*T⁴ + """ + R = 8.314462618 # kJ/kmol-K (universal gas constant) + a1, a2, a3, a4, a5, a6, a7 = coeffs + cp_over_R = a1 + a2*T + a3*T**2 + a4*T**3 + a5*T**4 + return cp_over_R * R + +def calculate_enthalpy_nasa(T, coeffs): + """ + Calculate enthalpy H(T) in kJ/kmol using NASA polynomial coefficients + H/RT = a1 + a2/2*T + a3/3*T² + a4/4*T³ + a5/5*T⁴ + a6/T + """ + R = 8.314462618 # kJ/kmol-K + a1, a2, a3, a4, a5, a6, a7 = coeffs + H_over_RT = a1 + a2/2*T + a3/3*T**2 + a4/4*T**3 + a5/5*T**4 + a6/T + return H_over_RT * R * T + +def calculate_adiabatic_flame_temperature(phi, T_in, c, tolerance=1e-3, max_iterations=100): + """ + Calculate adiabatic flame temperature for methane-air mixture with H2 addition + + Parameters: + phi: equivalence ratio + T_in: inlet temperature (K) + c: progress variable (0 to 1) + h2_fraction: fraction of H2 in the fuel mixture (0 to 1) + tolerance: convergence tolerance + max_iterations: maximum number of iterations + """ + # Calculate composition based on progress variable c + # For fuel mixture: (1 - h2_fraction) CH4 + h2_fraction H2 + n_CH4 = c*(0-0.6)+0.6 # Linear interpolation + n_H2 = c*(0-0.4)+0.4 # Linear interpolation + + # Oxygen requirement: CH4 needs 2O2, H2 needs 0.5O2 + + n_O2_air =c*(1.433-2.8)+2.8 + n_N2_air = 10.533 # Nitrogen scales with oxygen + + # Products composition + # CO2 from CH4 combustion only + n_CO2 = c*(0.605-0)+0 + # H2O from both CH4 and H2 combustion + n_H2O = c*(1.612-0)+0 + n_O2 = max(0, n_O2_air - (2 * n_CH4 + 0.5 * n_H2)) # Excess oxygen + n_N2 = n_N2_air + + # Enthalpy of formation (kJ/kmol) at 298.15K + hf_CH4 = -74850 + hf_H2 = 0 # H2 enthalpy of formation is 0 (element) + hf_O2 = 0 + hf_N2 = 0 + hf_CO2 = -393520 + hf_H2O = -241820 + + # Reference temperature + T_ref = 298.15 # K + + # Initial guess for adiabatic flame temperature + T_guess = 1800 # K + R = 8.314462618 # kJ/kmol-K + + for iteration in range(max_iterations): + # Get coefficients for current temperature guess (products) + coeffs_CO2 = get_nasa_coefficients('CO2', T_guess) + coeffs_H2O = get_nasa_coefficients('H2O', T_guess) + coeffs_O2_prod = get_nasa_coefficients('O2', T_guess) + coeffs_N2_prod = get_nasa_coefficients('N2', T_guess) + + # Get coefficients for inlet temperature (reactants) + coeffs_CH4_in = get_nasa_coefficients('CH4', T_in) + coeffs_H2_in = get_nasa_coefficients('H2', T_in) + coeffs_O2_in = get_nasa_coefficients('O2', T_in) + coeffs_N2_in = get_nasa_coefficients('N2', T_in) + + # Get coefficients for reference temperature + coeffs_CH4_ref = get_nasa_coefficients('CH4', T_ref) + coeffs_H2_ref = get_nasa_coefficients('H2', T_ref) + coeffs_O2_ref = get_nasa_coefficients('O2', T_ref) + coeffs_N2_ref = get_nasa_coefficients('N2', T_ref) + coeffs_CO2_ref = get_nasa_coefficients('CO2', T_ref) + coeffs_H2O_ref = get_nasa_coefficients('H2O', T_ref) + + # 1. Calculate Enthalpy of Reactants at inlet temperature + H_reactants = 0 + H_reactants += n_CH4 * (hf_CH4 + (calculate_enthalpy_nasa(T_in, coeffs_CH4_in) - calculate_enthalpy_nasa(T_ref, coeffs_CH4_ref))) + H_reactants += n_H2 * (hf_H2 + (calculate_enthalpy_nasa(T_in, coeffs_H2_in) - calculate_enthalpy_nasa(T_ref, coeffs_H2_ref))) + H_reactants += n_O2_air * (hf_O2 + (calculate_enthalpy_nasa(T_in, coeffs_O2_in) - calculate_enthalpy_nasa(T_ref, coeffs_O2_ref))) + H_reactants += n_N2_air * (hf_N2 + (calculate_enthalpy_nasa(T_in, coeffs_N2_in) - calculate_enthalpy_nasa(T_ref, coeffs_N2_ref))) + + # 2. Calculate Enthalpy of Products at guessed temperature + H_products = 0 + H_products += n_CO2 * (hf_CO2 + (calculate_enthalpy_nasa(T_guess, coeffs_CO2) - calculate_enthalpy_nasa(T_ref, coeffs_CO2_ref))) + H_products += n_H2O * (hf_H2O + (calculate_enthalpy_nasa(T_guess, coeffs_H2O) - calculate_enthalpy_nasa(T_ref, coeffs_H2O_ref))) + H_products += n_O2 * (hf_O2 + (calculate_enthalpy_nasa(T_guess, coeffs_O2_prod) - calculate_enthalpy_nasa(T_ref, coeffs_O2_ref))) + H_products += n_N2 * (hf_N2 + (calculate_enthalpy_nasa(T_guess, coeffs_N2_prod) - calculate_enthalpy_nasa(T_ref, coeffs_N2_ref))) + + # 3. Energy Balance Residual: H_reactants - H_products = 0 + residual = H_reactants - H_products + + # 4. Calculate total heat capacity of the product mixture at T_guess + cp_mix = (n_CO2 * calculate_cp_nasa(T_guess, coeffs_CO2) + + n_H2O * calculate_cp_nasa(T_guess, coeffs_H2O) + + n_O2 * calculate_cp_nasa(T_guess, coeffs_O2_prod) + + n_N2 * calculate_cp_nasa(T_guess, coeffs_N2_prod)) + + # 5. Newton-Raphson update: T_new = T_guess + residual / cp_mix + T_new = T_guess + residual / cp_mix + + # 6. Check for convergence + error = abs(T_new - T_guess) + + if error < tolerance: + return T_new + + T_guess = T_new + + return T_guess + +# P = rho*R*T +# 5 = 2.55 * R * 673 +# R = 0.0029 + +# ################################################################## # +# create a function for the initial progress variable # +# ################################################################## # +def initC(coord): + x = coord[0] + y = coord[1] + #z = coord[2] + #print("x,y = ",x," ",y) + C = 0.0 + # location where the flame should be + flame_x = 0.012 + if (x < flame_x): + C = 0.0 + else: + C = 1.0 + + return C + +# ################################################################## # +# loop over all vertices and set the species progress variable # +# ################################################################## # +def SetInitialSpecies(SU2Driver): + allCoords = SU2Driver.Coordinates() + iSPECIESSOLVER = SU2Driver.GetSolverIndices()['SPECIES'] + for iPoint in range(SU2Driver.GetNumberNodes() - SU2Driver.GetNumberHaloNodes()): + coord = allCoords.Get(iPoint) + C = initC(coord) + # now update the initial condition for the species + SU2Driver.Solution(iSPECIESSOLVER).Set(iPoint,0,C) + +def update_temperature(SU2Driver, iPoint): + # first, get the progress variable + iSPECIESSOLVER = SU2Driver.GetSolverIndices()['SPECIES'] + # Note: returns a list + C = SU2Driver.Solution(iSPECIESSOLVER)(iPoint,0) + h= SU2Driver.Solution(iSPECIESSOLVER)(iPoint,1) + #print("enthalpy",h) + + Tad=calculate_adiabatic_flame_temperature(phi,Tu,C) + #print("Adiabatic temp=",Tad) + y=-86.48*(C**4)+321.51*(C**3)-706.92*(C**2)+548.63*C-38.021 + Tad_f=Tad - y + T=(h-304721)/(1161) + Tad_f + + iFLOWSOLVER = SU2Driver.GetSolverIndices()['INC.FLOW'] + # the list with names + solindex = getsolvar(SU2Driver) + iTEMP = solindex.get("TEMPERATURE") + #print("Adiabatic temp=",T) + SU2Driver.Solution(iFLOWSOLVER).Set(iPoint,iTEMP,T) + + + +def zimont(SU2Driver, iPoint): + + iSSTSOLVER = SU2Driver.GetSolverIndices()['SST'] + tke, dissipation = SU2Driver.Solution(iSSTSOLVER)(iPoint) + + iSPECIESSOLVER = SU2Driver.GetSolverIndices()['SPECIES'] + # get the gradient of species_0 + gradc = SU2Driver.Gradient(iSPECIESSOLVER)(iPoint,0) + primindex = SU2Driver.GetPrimitiveIndices() + iDENSITY = primindex.get("DENSITY") + iMU = primindex.get("LAMINAR_VISCOSITY") + + h= SU2Driver.Solution(iSPECIESSOLVER)(iPoint,1) + h_scaled=h//1000000 + # laminar burning velocity of methane-air at phi=0.5, P=5 + + #Slu = 0.33745 + + rho = SU2Driver.Primitives()(iPoint,iDENSITY) + mu = SU2Driver.Primitives()(iPoint,iMU) + nu=mu/rho + # Turbulent Flamespeed Closure with Dinkelacker correction + up = np.sqrt((2.0/3.0) * tke ) + lt = (0.09**0.75) * (tke**1.5) / dissipation + Re = up*lt/nu + Le = 0.498 + if h_scaled>=-0.3: + Slu= -1.6508*(h_scaled**5)+1.3050*(h_scaled**4)+1.7592*(h_scaled**3)+1.1562*(h_scaled**2)+0.3947*(h_scaled)+0.0545 + Ut = Slu * (1.0 + (0.46/Le)*np.power(Re,0.25)*np.power(up/Slu,0.3)*np.power(Pu,0.2)) + else : + Ut = 0 + + norm_gradc = np.sqrt(gradc[0]*gradc[0] + gradc[1]*gradc[1]) + Sc = rho_u * Ut * norm_gradc + + return Sc + +def getsolvar(SU2Driver): + primindex = SU2Driver.GetPrimitiveIndices() + iFLOWSOLVER = SU2Driver.GetSolverIndices()['INC.FLOW'] + nVars = SU2Driver.Solution(iFLOWSOLVER).Shape()[1] + varindex = primindex.copy() + for prim in varindex.copy(): + if varindex[prim] >=nVars: + del varindex[prim] + varindex = dict(sorted(varindex.items(), key=lambda item: item[1])) + return varindex + +def ApplyScalar(driver,marker_ids): + for marker_id in marker_ids: + if marker_id < 0: + continue + + #h1=[] + #h2=[] + #print('Marker id',marker_id) + for i_vertex in range(driver.GetNumberMarkerNodes(marker_id)): + iSPECIESSOLVER = driver.GetSolverIndices()['SPECIES'] + C=driver.Solution(iSPECIESSOLVER)(i_vertex,0) + Tad=calculate_adiabatic_flame_temperature(phi,Tu,C) + y=-86.48*(C**4)+321.51*(C**3)-706.92*(C**2)+548.63*C-38.021 #Correction factor + Tad_f=Tad - y + h1=1 + h2=304721+1161*(Twall-Tad_f) + hf=[h1,h2] + driver.SetMarkerCustomScalar(marker_id, i_vertex, hf) + +def main(): + """ + Run the flow solver with a custom inlet (function of time and space). + """ + comm = MPI.COMM_WORLD + #comm = 0 + rank = comm.Get_rank() + # Initialize the primal driver of SU2, this includes solver preprocessing. + try: + driver = pysu2.CSinglezoneDriver('psi.cfg', 1, comm) + except TypeError as exception: + print('A TypeError occured in pysu2.CSinglezoneDriver : ', exception) + raise + + print("\n------------------------------ Begin Solver -----------------------------") + sys.stdout.flush() + + nDim = driver.GetNumberDimensions() + + # index to the flow solver + # C.FLOW + # INC.FLOW + # HEAT + # FLAMELET + # SPECIES + # SA + # SST + iFLOWSOLVER = driver.GetSolverIndices()['INC.FLOW'] + iSPECIESSOLVER = driver.GetSolverIndices()['SPECIES'] + iSSTSOLVER = driver.GetSolverIndices()['SST'] + # all the indices and the map to the names of the primitives + primindex = driver.GetPrimitiveIndices() + nElem = driver.GetNumberElements() + nVars = driver.Solution(iFLOWSOLVER).Shape()[1] + nVarsSpecies = driver.Solution(iSPECIESSOLVER).Shape()[1] + nVarsTurb = driver.Solution(iSSTSOLVER).Shape()[1] + + solindex = getsolvar(driver) + iTEMP = solindex.get("TEMPERATURE") + + if rank == 0: + print("Dimensions of the problem = ",nDim) + print("index of flow solver = ",iFLOWSOLVER) + print("index of turbulence solver = ",iSSTSOLVER) + print("indices of primitives=",primindex) + print("number of primitives:",len(primindex)) + print("number of elements:",nElem) + print("number of flow solver variables:",nVars) + print("number of species solver variables:",nVarsSpecies) + print("number of turbulence solver variables:",nVarsTurb) + sys.stdout.flush() + + + with open('psi.cfg') as f: + if 'RESTART_SOL= YES' in f.read(): + if rank == 0: + print("restarting from file") + else: + # We can set an initial condition by calling this function: + if rank == 0: + print("Using user defined initial condition.") + SetInitialSpecies(driver) + + sys.stdout.flush() + + + all_marker_ids = driver.GetMarkerIndices() + marker_names = ['wall_top'] + #print('Marker id trial',all_marker_ids) + marker_ids = [] + for name in marker_names: + marker_ids.append(all_marker_ids[name] if name in all_marker_ids else -1) + + # run 500 iterations + for inner_iter in range(2000): + + if (rank==0): + print("python iteration ", inner_iter) + #ApplyScalar(driver,marker_ids) + driver.Preprocess(inner_iter) + driver.Run() + + Source = driver.UserDefinedSource(iSPECIESSOLVER) + + # set the source term, per point + for i_node in range(driver.GetNumberNodes() - driver.GetNumberHaloNodes()): + # add source term: + # default TFC of Zimont: rho*Sc = rho_u * U_t * grad(c) + S = zimont(driver, i_node) + Source.Set(i_node, 0, S) + + # for the update of temperature, we need to update also the halo nodes + #for i_node in range(driver.GetNumberNodes()): + # # set the temperature to T = c*Tf + (1-c)*Tu + update_temperature(driver, i_node) + + driver.Postprocess() + driver.Update() + # Monitor the solver and output solution to file if required. + #driver.Monitor(inner_iter) + # Output the solution to file + driver.Output(inner_iter) + + # Finalize the solver and exit cleanly. + driver.Finalize() + +if __name__ == '__main__': + main() diff --git a/TestCases/py_wrapper/psi_heatloss/vol_solution/zone_0/symmetry.vtu-934412289-428597.lock b/TestCases/py_wrapper/psi_heatloss/vol_solution/zone_0/symmetry.vtu-934412289-428597.lock new file mode 100644 index 0000000000000000000000000000000000000000..1b1cb4d44c57c2d7a5122870fa6ac3e62ff7e94e GIT binary patch literal 8 KcmZQzfB*mh2mk>9 literal 0 HcmV?d00001 From 39ad1986fb459c6b4d9107a53dceb9100b1d5319 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Thu, 4 Dec 2025 21:42:34 +0100 Subject: [PATCH 11/31] run precommit --- Common/include/option_structure.inl | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/Common/include/option_structure.inl b/Common/include/option_structure.inl index d3592564360e..2f149d3b97c6 100644 --- a/Common/include/option_structure.inl +++ b/Common/include/option_structure.inl @@ -1317,7 +1317,10 @@ class COptionWallSpecies : public COptionBase { COptionWallSpecies(string option_field_name, unsigned short& nMarker_Wall_Species, string*& Marker_Wall_Species, unsigned short**& option_field, const map m, su2double**& value, unsigned short& nSpecies_per_Wall) - : size(nMarker_Wall_Species), marker(Marker_Wall_Species), field(option_field), value(value), + : size(nMarker_Wall_Species), + marker(Marker_Wall_Species), + field(option_field), + value(value), nSpecies_per_Wall(nSpecies_per_Wall) { this->name = option_field_name; this->m = m; From ff27515d44f897d49122cf434904be04be6a1735 Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 17 Dec 2025 10:03:47 +0100 Subject: [PATCH 12/31] Potential fix for code scanning alert no. 5867: Unused local variable Co-authored-by: Copilot Autofix powered by AI <62310815+github-advanced-security[bot]@users.noreply.github.com> --- TestCases/py_wrapper/psi_heatloss/run.py | 1 - 1 file changed, 1 deletion(-) diff --git a/TestCases/py_wrapper/psi_heatloss/run.py b/TestCases/py_wrapper/psi_heatloss/run.py index 92511a31a99a..80c6eb95a209 100755 --- a/TestCases/py_wrapper/psi_heatloss/run.py +++ b/TestCases/py_wrapper/psi_heatloss/run.py @@ -151,7 +151,6 @@ def calculate_adiabatic_flame_temperature(phi, T_in, c, tolerance=1e-3, max_iter # Initial guess for adiabatic flame temperature T_guess = 1800 # K - R = 8.314462618 # kJ/kmol-K for iteration in range(max_iterations): # Get coefficients for current temperature guess (products) From 40d5b99fee09be160a3a04ef4a9e31b324f5d071 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Wed, 17 Dec 2025 11:06:57 +0100 Subject: [PATCH 13/31] cleanup --- TestCases/py_wrapper/psi_heatloss/psi.cfg | 176 ------- TestCases/py_wrapper/psi_heatloss/run.py | 443 ------------------ .../zone_0/symmetry.vtu-934412289-428597.lock | Bin 8 -> 0 bytes .../species2_primitiveVenturi_flux_value.cfg | 131 ------ .../species3_primitiveVenturi_flux_value.cfg | 21 +- config_template.cfg | 1 + 6 files changed, 12 insertions(+), 760 deletions(-) delete mode 100755 TestCases/py_wrapper/psi_heatloss/psi.cfg delete mode 100755 TestCases/py_wrapper/psi_heatloss/run.py delete mode 100644 TestCases/py_wrapper/psi_heatloss/vol_solution/zone_0/symmetry.vtu-934412289-428597.lock delete mode 100644 TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_flux_value.cfg diff --git a/TestCases/py_wrapper/psi_heatloss/psi.cfg b/TestCases/py_wrapper/psi_heatloss/psi.cfg deleted file mode 100755 index 6d1cb64275ed..000000000000 --- a/TestCases/py_wrapper/psi_heatloss/psi.cfg +++ /dev/null @@ -1,176 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% % -% SU2 configuration file % -% Case description: Species mixing with 3 species, i.e. 2 transport equations % -% Author: T. Kattmann % -% Institution: Bosch Thermotechniek B.V. % -% Date: 2021/10/14 % -% File Version 7.5.1 "Blackbird" % -% % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% -% -SOLVER= INC_RANS -KIND_TURB_MODEL= SST -SST_OPTIONS= V2003m - -RESTART_SOL= YES - -MGLEVEL= 0 -% - -% temperature and density -FREESTREAM_TEMPERATURE = 673 -FREESTREAM_DENSITY = 2.50 -% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% -% -%INC_DENSITY_MODEL= CONSTANT -INC_DENSITY_MODEL= VARIABLE -INC_DENSITY_INIT= 2.50 -% -INC_VELOCITY_INIT= (40.00, 0.0, 0.0 ) -% -INC_ENERGY_EQUATION= YES -INC_TEMPERATURE_INIT= 673.0 -% -INC_NONDIM= DIMENSIONAL -% -% -------------------- FLUID PROPERTIES ------------------------------------- % -% -%FLUID_MODEL= CONSTANT_DENSITY -FLUID_MODEL= INC_IDEAL_GAS -% -CONDUCTIVITY_MODEL= CONSTANT_CONDUCTIVITY -THERMAL_CONDUCTIVITY_CONSTANT= 0.0357 -% -PRANDTL_LAM= 0.72 -TURBULENT_CONDUCTIVITY_MODEL= NONE -PRANDTL_TURB= 0.90 -% -VISCOSITY_MODEL= SUTHERLAND -MU_CONSTANT= 1.716E-5 -MU_REF = 1.716e-5 -MU_T_REF= 273.15 -SUTHERLAND_CONSTANT = 110.4 - -SPECIFIC_HEAT_CP = 1150 -% -% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% -% -MARKER_HEATFLUX= ( wall_pipe,0.0, wall_interior, 0.0,wall_top, 0,wall_side,0 ) - -% note, case is axisymmetric -MARKER_SYM= ( symmetry ) -AXISYMMETRIC= YES -SPECIFIED_INLET_PROFILE= NO -INLET_FILENAME= inlet.dat -INC_INLET_TYPE= VELOCITY_INLET -INC_INLET_DAMPING= 0.01 -MARKER_INLET= ( inlet, 673, 40.0, 1.0, 0.0, 0.0) -MARKER_INLET_TURBULENT = (inlet, 0.10, 15) -MARKER_INLET_SPECIES= (inlet, 0.0, 304721) -INLET_MATCHING_TOLERANCE = 1e-3 -MARKER_PYTHON_CUSTOM= (wall_top) - -MARKER_WALL_SPECIES= wall_side, FLUX,0,VALUE,304721, wall_pipe, FLUX,0 , VALUE,304721 ,wall_interior , FLUX,0 , VALUE,304721 ,wall_top, FLUX,0,VALUE,-250000 - -MARKER_SPECIES_STRONG_BC=wall_top,wall_side,wall_pipe,wall_interior -% -INC_OUTLET_TYPE= PRESSURE_OUTLET -INC_OUTLET_DAMPING= 0.01 -MARKER_OUTLET= ( outlet, 0.0 ) -% -% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% -% -NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES -% -%CFL_NUMBER= 25.0 -CFL_NUMBER= 10.0 -CFL_REDUCTION_SPECIES= 1.0 -CFL_REDUCTION_TURB= 1.0 -CFL_ADAPT= NO -CFL_ADAPT_PARAM= ( 0.95, 1.01, 1.0, 250, 1.0e-4, 0) -% -% Run commented Iter for good results -ITER= 1 -% -% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% -% -LINEAR_SOLVER= FGMRES -LINEAR_SOLVER_PREC= ILU -LINEAR_SOLVER_ERROR= 1E-10 -LINEAR_SOLVER_ITER= 10 - -% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% -% -CONV_NUM_METHOD_FLOW= FDS -MUSCL_FLOW= NO -SLOPE_LIMITER_FLOW = NONE -TIME_DISCRE_FLOW= EULER_IMPLICIT -% -% -------------------- SCALAR TRANSPORT ---------------------------------------% -% -KIND_SCALAR_MODEL= SPECIES_TRANSPORT -%DIFFUSIVITY_MODEL= CONSTANT_SCHMIDT -DIFFUSIVITY_MODEL= CONSTANT_DIFFUSIVITY -SCHMIDT_NUMBER_LAMINAR= 0.5 -DIFFUSIVITY_CONSTANT= 0.000148 - -% according to the paper -SCHMIDT_NUMBER_TURBULENT= 0.7 -% -CONV_NUM_METHOD_SPECIES= BOUNDED_SCALAR -MUSCL_SPECIES= NO -SLOPE_LIMITER_SPECIES = NONE -% -TIME_DISCRE_SPECIES= EULER_IMPLICIT -% -SPECIES_INIT= 0.0, 200000 - -SPECIES_CLIPPING= NO -SPECIES_CLIPPING_MIN= 0.0, -1.5e+6 -SPECIES_CLIPPING_MAX= 1.0, 304721 - -% -% -------------------- TURBULENT TRANSPORT ---------------------------------------% -% -CONV_NUM_METHOD_TURB= BOUNDED_SCALAR -MUSCL_TURB= NO - -% -% --------------------------- CONVERGENCE PARAMETERS --------------------------% -% -CONV_FIELD= RMS_PRESSURE, RMS_VELOCITY-X, RMS_VELOCITY-Y, RMS_TKE, RMS_SPECIES -CONV_RESIDUAL_MINVAL= -18 -CONV_STARTITER= 10 -% -% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% -% -MESH_FILENAME= psi.su2 -%MESH_FILENAME= psi_medium.su2 -%MESH_FILENAME= psi_fine.su2 -% -SCREEN_OUTPUT= INNER_ITER WALL_TIME \ - RMS_PRESSURE RMS_VELOCITY-X RMS_VELOCITY-Y RMS_TKE RMS_DISSIPATION RMS_SPECIES_0 RMS_SPECIES_1 \ - LINSOL_ITER LINSOL_RESIDUAL \ - LINSOL_ITER_TURB LINSOL_RESIDUAL_TURB \ - LINSOL_ITER_SPECIES LINSOL_RESIDUAL_SPECIES \ - SURFACE_SPECIES_0 -SCREEN_WRT_FREQ_INNER= 1 -% -HISTORY_OUTPUT= ITER RMS_RES LINSOL SPECIES_COEFF SPECIES_COEFF_SURF -CONV_FILENAME= history -MARKER_ANALYZE= gas_inlet, air_axial_inlet, outlet -MARKER_ANALYZE_AVERAGE= AREA -% -OUTPUT_FILES= RESTART_ASCII, PARAVIEW_MULTIBLOCK -VOLUME_OUTPUT= RESIDUAL, PRIMITIVE, RANK, SPECIES_UDS_0 -OUTPUT_WRT_FREQ= 100 -% -READ_BINARY_RESTART= NO -RESTART_FILENAME= restart -SOLUTION_FILENAME= solution -% -WRT_PERFORMANCE= YES -PYTHON_CUSTOM_SOURCE= YES diff --git a/TestCases/py_wrapper/psi_heatloss/run.py b/TestCases/py_wrapper/psi_heatloss/run.py deleted file mode 100755 index 80c6eb95a209..000000000000 --- a/TestCases/py_wrapper/psi_heatloss/run.py +++ /dev/null @@ -1,443 +0,0 @@ -#!/usr/bin/env python - -## \file run.py -# \brief turbulent premixed dump combustor simulation (PSI flame) -# phi=0.5, methane-air, U=40 m/s -# \version 8.1.0 "Harrier" -# -# SU2 Project Website: https://su2code.github.io -# -# The SU2 Project is maintained by the SU2 Foundation -# (http://su2foundation.org) -# -# Copyright 2012-2024, 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 . - -import sys -import pysu2 -import numpy as np -from mpi4py import MPI - -# unburnt temperature of the propane-air mixture -# flame temperature of the methane-air mixture (phi=0.5, P=5) -Tf = 1803 -Twall=1200 -Tu = 673.0 -Pu = 5.0 -phi = 0.5 -# unburnt density at P=5 -rho_u = 2.50 -# unburnt thermal conductivity of methane-air at phi=0.5 (P=5) -k_u = 0.0523 -# unburnt heat capacity of methane-air at phi=0.5 (P=1) -cp_u = 1311.0 - -NASA_COEFFICIENTS = { - 'CH4': { - 'low_coeffs': [5.14911468, -0.013662201, 4.91454E-05, -4.84247E-08, 1.66603E-11, -10246.5983, -4.63848842], - 'high_coeffs': [1.65326226, 0.01002631, -3.31661E-06, 5.36483E-10, -3.14697E-14, -10009.5936, 9.90506283] - }, - 'O2': { - 'low_coeffs': [3.78245636, -0.002996734, 9.8473E-06, -9.6813E-09, 3.24373E-12, -1063.94356, 3.65767573], - 'high_coeffs': [3.66096083, 0.000656366, -1.41149E-07, 2.05798E-11, -1.29913E-15, -1215.97725, 3.41536184] - }, - 'N2': { - 'low_coeffs': [3.53100528, -0.000123661, -5.02999E-07, 2.43531E-09, -1.40881E-12, -1046.97628, 2.96747038], - 'high_coeffs': [2.95257637, 0.0013969, -4.92632E-07, 7.8601E-11, -4.60755E-15, -923.948688, 5.87188762] - }, - 'CO2': { - 'low_coeffs': [2.356813, 0.00898413, -7.12206E-06, 2.4573E-09, -1.42885E-13, -48371.971, 9.9009035], - 'high_coeffs': [4.6365111, 0.002741457, -9.95898E-07, 1.60387E-10, -9.16199E-15, -49024.904, -1.9348955] - }, - 'H2O': { - 'low_coeffs': [4.1986352, -0.002036402, 6.52034E-06, -5.48793E-09, 1.77197E-12, -30293.726, -0.84900901], - 'high_coeffs': [2.6770389, 0.002973182, -7.73769E-07, 9.44335E-11, -4.269E-15, -29885.894, 6.88255] - }, - 'H2': { - 'low_coeffs': [2.34433112, 0.007980521, -1.94782E-05, 2.01572E-08, -7.37612E-12, -917.935173, 0.683010238], - 'high_coeffs': [2.93286575, 0.000826608, -1.46402E-07, 1.541E-11, -6.88805E-16, -813.065581, -1.02432865] - } -} - -# Temperature range boundaries -T_LOW = 300.0 # K -T_HIGH = 1000.0 # K - -def get_nasa_coefficients(species, T): - """ - Return the appropriate coefficients for the given temperature - """ - if species not in NASA_COEFFICIENTS: - raise ValueError(f"Species {species} not found in database") - - if T <= T_HIGH: - return NASA_COEFFICIENTS[species]['low_coeffs'] - else: - return NASA_COEFFICIENTS[species]['high_coeffs'] - - -def calculate_cp_nasa(T, coeffs): - """ - Calculate Cp in kJ/kmol-K using NASA polynomial coefficients - Cp/R = a1 + a2*T + a3*T² + a4*T³ + a5*T⁴ - """ - R = 8.314462618 # kJ/kmol-K (universal gas constant) - a1, a2, a3, a4, a5, a6, a7 = coeffs - cp_over_R = a1 + a2*T + a3*T**2 + a4*T**3 + a5*T**4 - return cp_over_R * R - -def calculate_enthalpy_nasa(T, coeffs): - """ - Calculate enthalpy H(T) in kJ/kmol using NASA polynomial coefficients - H/RT = a1 + a2/2*T + a3/3*T² + a4/4*T³ + a5/5*T⁴ + a6/T - """ - R = 8.314462618 # kJ/kmol-K - a1, a2, a3, a4, a5, a6, a7 = coeffs - H_over_RT = a1 + a2/2*T + a3/3*T**2 + a4/4*T**3 + a5/5*T**4 + a6/T - return H_over_RT * R * T - -def calculate_adiabatic_flame_temperature(phi, T_in, c, tolerance=1e-3, max_iterations=100): - """ - Calculate adiabatic flame temperature for methane-air mixture with H2 addition - - Parameters: - phi: equivalence ratio - T_in: inlet temperature (K) - c: progress variable (0 to 1) - h2_fraction: fraction of H2 in the fuel mixture (0 to 1) - tolerance: convergence tolerance - max_iterations: maximum number of iterations - """ - # Calculate composition based on progress variable c - # For fuel mixture: (1 - h2_fraction) CH4 + h2_fraction H2 - n_CH4 = c*(0-0.6)+0.6 # Linear interpolation - n_H2 = c*(0-0.4)+0.4 # Linear interpolation - - # Oxygen requirement: CH4 needs 2O2, H2 needs 0.5O2 - - n_O2_air =c*(1.433-2.8)+2.8 - n_N2_air = 10.533 # Nitrogen scales with oxygen - - # Products composition - # CO2 from CH4 combustion only - n_CO2 = c*(0.605-0)+0 - # H2O from both CH4 and H2 combustion - n_H2O = c*(1.612-0)+0 - n_O2 = max(0, n_O2_air - (2 * n_CH4 + 0.5 * n_H2)) # Excess oxygen - n_N2 = n_N2_air - - # Enthalpy of formation (kJ/kmol) at 298.15K - hf_CH4 = -74850 - hf_H2 = 0 # H2 enthalpy of formation is 0 (element) - hf_O2 = 0 - hf_N2 = 0 - hf_CO2 = -393520 - hf_H2O = -241820 - - # Reference temperature - T_ref = 298.15 # K - - # Initial guess for adiabatic flame temperature - T_guess = 1800 # K - - for iteration in range(max_iterations): - # Get coefficients for current temperature guess (products) - coeffs_CO2 = get_nasa_coefficients('CO2', T_guess) - coeffs_H2O = get_nasa_coefficients('H2O', T_guess) - coeffs_O2_prod = get_nasa_coefficients('O2', T_guess) - coeffs_N2_prod = get_nasa_coefficients('N2', T_guess) - - # Get coefficients for inlet temperature (reactants) - coeffs_CH4_in = get_nasa_coefficients('CH4', T_in) - coeffs_H2_in = get_nasa_coefficients('H2', T_in) - coeffs_O2_in = get_nasa_coefficients('O2', T_in) - coeffs_N2_in = get_nasa_coefficients('N2', T_in) - - # Get coefficients for reference temperature - coeffs_CH4_ref = get_nasa_coefficients('CH4', T_ref) - coeffs_H2_ref = get_nasa_coefficients('H2', T_ref) - coeffs_O2_ref = get_nasa_coefficients('O2', T_ref) - coeffs_N2_ref = get_nasa_coefficients('N2', T_ref) - coeffs_CO2_ref = get_nasa_coefficients('CO2', T_ref) - coeffs_H2O_ref = get_nasa_coefficients('H2O', T_ref) - - # 1. Calculate Enthalpy of Reactants at inlet temperature - H_reactants = 0 - H_reactants += n_CH4 * (hf_CH4 + (calculate_enthalpy_nasa(T_in, coeffs_CH4_in) - calculate_enthalpy_nasa(T_ref, coeffs_CH4_ref))) - H_reactants += n_H2 * (hf_H2 + (calculate_enthalpy_nasa(T_in, coeffs_H2_in) - calculate_enthalpy_nasa(T_ref, coeffs_H2_ref))) - H_reactants += n_O2_air * (hf_O2 + (calculate_enthalpy_nasa(T_in, coeffs_O2_in) - calculate_enthalpy_nasa(T_ref, coeffs_O2_ref))) - H_reactants += n_N2_air * (hf_N2 + (calculate_enthalpy_nasa(T_in, coeffs_N2_in) - calculate_enthalpy_nasa(T_ref, coeffs_N2_ref))) - - # 2. Calculate Enthalpy of Products at guessed temperature - H_products = 0 - H_products += n_CO2 * (hf_CO2 + (calculate_enthalpy_nasa(T_guess, coeffs_CO2) - calculate_enthalpy_nasa(T_ref, coeffs_CO2_ref))) - H_products += n_H2O * (hf_H2O + (calculate_enthalpy_nasa(T_guess, coeffs_H2O) - calculate_enthalpy_nasa(T_ref, coeffs_H2O_ref))) - H_products += n_O2 * (hf_O2 + (calculate_enthalpy_nasa(T_guess, coeffs_O2_prod) - calculate_enthalpy_nasa(T_ref, coeffs_O2_ref))) - H_products += n_N2 * (hf_N2 + (calculate_enthalpy_nasa(T_guess, coeffs_N2_prod) - calculate_enthalpy_nasa(T_ref, coeffs_N2_ref))) - - # 3. Energy Balance Residual: H_reactants - H_products = 0 - residual = H_reactants - H_products - - # 4. Calculate total heat capacity of the product mixture at T_guess - cp_mix = (n_CO2 * calculate_cp_nasa(T_guess, coeffs_CO2) + - n_H2O * calculate_cp_nasa(T_guess, coeffs_H2O) + - n_O2 * calculate_cp_nasa(T_guess, coeffs_O2_prod) + - n_N2 * calculate_cp_nasa(T_guess, coeffs_N2_prod)) - - # 5. Newton-Raphson update: T_new = T_guess + residual / cp_mix - T_new = T_guess + residual / cp_mix - - # 6. Check for convergence - error = abs(T_new - T_guess) - - if error < tolerance: - return T_new - - T_guess = T_new - - return T_guess - -# P = rho*R*T -# 5 = 2.55 * R * 673 -# R = 0.0029 - -# ################################################################## # -# create a function for the initial progress variable # -# ################################################################## # -def initC(coord): - x = coord[0] - y = coord[1] - #z = coord[2] - #print("x,y = ",x," ",y) - C = 0.0 - # location where the flame should be - flame_x = 0.012 - if (x < flame_x): - C = 0.0 - else: - C = 1.0 - - return C - -# ################################################################## # -# loop over all vertices and set the species progress variable # -# ################################################################## # -def SetInitialSpecies(SU2Driver): - allCoords = SU2Driver.Coordinates() - iSPECIESSOLVER = SU2Driver.GetSolverIndices()['SPECIES'] - for iPoint in range(SU2Driver.GetNumberNodes() - SU2Driver.GetNumberHaloNodes()): - coord = allCoords.Get(iPoint) - C = initC(coord) - # now update the initial condition for the species - SU2Driver.Solution(iSPECIESSOLVER).Set(iPoint,0,C) - -def update_temperature(SU2Driver, iPoint): - # first, get the progress variable - iSPECIESSOLVER = SU2Driver.GetSolverIndices()['SPECIES'] - # Note: returns a list - C = SU2Driver.Solution(iSPECIESSOLVER)(iPoint,0) - h= SU2Driver.Solution(iSPECIESSOLVER)(iPoint,1) - #print("enthalpy",h) - - Tad=calculate_adiabatic_flame_temperature(phi,Tu,C) - #print("Adiabatic temp=",Tad) - y=-86.48*(C**4)+321.51*(C**3)-706.92*(C**2)+548.63*C-38.021 - Tad_f=Tad - y - T=(h-304721)/(1161) + Tad_f - - iFLOWSOLVER = SU2Driver.GetSolverIndices()['INC.FLOW'] - # the list with names - solindex = getsolvar(SU2Driver) - iTEMP = solindex.get("TEMPERATURE") - #print("Adiabatic temp=",T) - SU2Driver.Solution(iFLOWSOLVER).Set(iPoint,iTEMP,T) - - - -def zimont(SU2Driver, iPoint): - - iSSTSOLVER = SU2Driver.GetSolverIndices()['SST'] - tke, dissipation = SU2Driver.Solution(iSSTSOLVER)(iPoint) - - iSPECIESSOLVER = SU2Driver.GetSolverIndices()['SPECIES'] - # get the gradient of species_0 - gradc = SU2Driver.Gradient(iSPECIESSOLVER)(iPoint,0) - primindex = SU2Driver.GetPrimitiveIndices() - iDENSITY = primindex.get("DENSITY") - iMU = primindex.get("LAMINAR_VISCOSITY") - - h= SU2Driver.Solution(iSPECIESSOLVER)(iPoint,1) - h_scaled=h//1000000 - # laminar burning velocity of methane-air at phi=0.5, P=5 - - #Slu = 0.33745 - - rho = SU2Driver.Primitives()(iPoint,iDENSITY) - mu = SU2Driver.Primitives()(iPoint,iMU) - nu=mu/rho - # Turbulent Flamespeed Closure with Dinkelacker correction - up = np.sqrt((2.0/3.0) * tke ) - lt = (0.09**0.75) * (tke**1.5) / dissipation - Re = up*lt/nu - Le = 0.498 - if h_scaled>=-0.3: - Slu= -1.6508*(h_scaled**5)+1.3050*(h_scaled**4)+1.7592*(h_scaled**3)+1.1562*(h_scaled**2)+0.3947*(h_scaled)+0.0545 - Ut = Slu * (1.0 + (0.46/Le)*np.power(Re,0.25)*np.power(up/Slu,0.3)*np.power(Pu,0.2)) - else : - Ut = 0 - - norm_gradc = np.sqrt(gradc[0]*gradc[0] + gradc[1]*gradc[1]) - Sc = rho_u * Ut * norm_gradc - - return Sc - -def getsolvar(SU2Driver): - primindex = SU2Driver.GetPrimitiveIndices() - iFLOWSOLVER = SU2Driver.GetSolverIndices()['INC.FLOW'] - nVars = SU2Driver.Solution(iFLOWSOLVER).Shape()[1] - varindex = primindex.copy() - for prim in varindex.copy(): - if varindex[prim] >=nVars: - del varindex[prim] - varindex = dict(sorted(varindex.items(), key=lambda item: item[1])) - return varindex - -def ApplyScalar(driver,marker_ids): - for marker_id in marker_ids: - if marker_id < 0: - continue - - #h1=[] - #h2=[] - #print('Marker id',marker_id) - for i_vertex in range(driver.GetNumberMarkerNodes(marker_id)): - iSPECIESSOLVER = driver.GetSolverIndices()['SPECIES'] - C=driver.Solution(iSPECIESSOLVER)(i_vertex,0) - Tad=calculate_adiabatic_flame_temperature(phi,Tu,C) - y=-86.48*(C**4)+321.51*(C**3)-706.92*(C**2)+548.63*C-38.021 #Correction factor - Tad_f=Tad - y - h1=1 - h2=304721+1161*(Twall-Tad_f) - hf=[h1,h2] - driver.SetMarkerCustomScalar(marker_id, i_vertex, hf) - -def main(): - """ - Run the flow solver with a custom inlet (function of time and space). - """ - comm = MPI.COMM_WORLD - #comm = 0 - rank = comm.Get_rank() - # Initialize the primal driver of SU2, this includes solver preprocessing. - try: - driver = pysu2.CSinglezoneDriver('psi.cfg', 1, comm) - except TypeError as exception: - print('A TypeError occured in pysu2.CSinglezoneDriver : ', exception) - raise - - print("\n------------------------------ Begin Solver -----------------------------") - sys.stdout.flush() - - nDim = driver.GetNumberDimensions() - - # index to the flow solver - # C.FLOW - # INC.FLOW - # HEAT - # FLAMELET - # SPECIES - # SA - # SST - iFLOWSOLVER = driver.GetSolverIndices()['INC.FLOW'] - iSPECIESSOLVER = driver.GetSolverIndices()['SPECIES'] - iSSTSOLVER = driver.GetSolverIndices()['SST'] - # all the indices and the map to the names of the primitives - primindex = driver.GetPrimitiveIndices() - nElem = driver.GetNumberElements() - nVars = driver.Solution(iFLOWSOLVER).Shape()[1] - nVarsSpecies = driver.Solution(iSPECIESSOLVER).Shape()[1] - nVarsTurb = driver.Solution(iSSTSOLVER).Shape()[1] - - solindex = getsolvar(driver) - iTEMP = solindex.get("TEMPERATURE") - - if rank == 0: - print("Dimensions of the problem = ",nDim) - print("index of flow solver = ",iFLOWSOLVER) - print("index of turbulence solver = ",iSSTSOLVER) - print("indices of primitives=",primindex) - print("number of primitives:",len(primindex)) - print("number of elements:",nElem) - print("number of flow solver variables:",nVars) - print("number of species solver variables:",nVarsSpecies) - print("number of turbulence solver variables:",nVarsTurb) - sys.stdout.flush() - - - with open('psi.cfg') as f: - if 'RESTART_SOL= YES' in f.read(): - if rank == 0: - print("restarting from file") - else: - # We can set an initial condition by calling this function: - if rank == 0: - print("Using user defined initial condition.") - SetInitialSpecies(driver) - - sys.stdout.flush() - - - all_marker_ids = driver.GetMarkerIndices() - marker_names = ['wall_top'] - #print('Marker id trial',all_marker_ids) - marker_ids = [] - for name in marker_names: - marker_ids.append(all_marker_ids[name] if name in all_marker_ids else -1) - - # run 500 iterations - for inner_iter in range(2000): - - if (rank==0): - print("python iteration ", inner_iter) - #ApplyScalar(driver,marker_ids) - driver.Preprocess(inner_iter) - driver.Run() - - Source = driver.UserDefinedSource(iSPECIESSOLVER) - - # set the source term, per point - for i_node in range(driver.GetNumberNodes() - driver.GetNumberHaloNodes()): - # add source term: - # default TFC of Zimont: rho*Sc = rho_u * U_t * grad(c) - S = zimont(driver, i_node) - Source.Set(i_node, 0, S) - - # for the update of temperature, we need to update also the halo nodes - #for i_node in range(driver.GetNumberNodes()): - # # set the temperature to T = c*Tf + (1-c)*Tu - update_temperature(driver, i_node) - - driver.Postprocess() - driver.Update() - # Monitor the solver and output solution to file if required. - #driver.Monitor(inner_iter) - # Output the solution to file - driver.Output(inner_iter) - - # Finalize the solver and exit cleanly. - driver.Finalize() - -if __name__ == '__main__': - main() diff --git a/TestCases/py_wrapper/psi_heatloss/vol_solution/zone_0/symmetry.vtu-934412289-428597.lock b/TestCases/py_wrapper/psi_heatloss/vol_solution/zone_0/symmetry.vtu-934412289-428597.lock deleted file mode 100644 index 1b1cb4d44c57c2d7a5122870fa6ac3e62ff7e94e..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 8 KcmZQzfB*mh2mk>9 diff --git a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_flux_value.cfg b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_flux_value.cfg deleted file mode 100644 index f601cefc4668..000000000000 --- a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_flux_value.cfg +++ /dev/null @@ -1,131 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% % -% SU2 configuration file % -% Case description: Species mixing with 2 species, i.e. 1 transport equations % -% Author: T. Kattmann % -% Institution: Bosch Thermotechniek B.V. % -% Date: 2021/10/14 % -% File Version 8.3.0 "Harrier" % -% % -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% -% -SOLVER= INC_RANS -KIND_TURB_MODEL= SST -% -% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% -% -INC_DENSITY_MODEL= CONSTANT -INC_DENSITY_INIT= 1.1766 -% -INC_VELOCITY_INIT= ( 1.00, 0.0, 0.0 ) -% -INC_ENERGY_EQUATION= YES -INC_TEMPERATURE_INIT= 300.0 -% -INC_NONDIM= DIMENSIONAL -% -% -------------------- FLUID PROPERTIES ------------------------------------- % -% -FLUID_MODEL= CONSTANT_DENSITY -% -CONDUCTIVITY_MODEL= CONSTANT_CONDUCTIVITY -THERMAL_CONDUCTIVITY_CONSTANT= 0.0357 -% -PRANDTL_LAM= 0.72 -TURBULENT_CONDUCTIVITY_MODEL= NONE -PRANDTL_TURB= 0.90 -% -VISCOSITY_MODEL= CONSTANT_VISCOSITY -MU_CONSTANT= 1.716E-5 -% -% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% -% -MARKER_HEATFLUX= ( wall, 0.0 ) -MARKER_SYM= ( axis ) -% -INC_INLET_TYPE= VELOCITY_INLET VELOCITY_INLET -MARKER_INLET= ( gas_inlet, 300, 1.0, 1.0, 0.0, 0.0,\ - air_axial_inlet, 300, 1.0, 0.0, -1.0, 0.0 ) -MARKER_INLET_SPECIES= (gas_inlet, 0.1,\ - air_axial_inlet, 0.0 ) -MARKER_WALL_SPECIES=(axis,FLUX,0.1,wall,VALUE,0.5) -% -INC_OUTLET_TYPE= PRESSURE_OUTLET -MARKER_OUTLET= ( outlet, 0.0) -% -% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% -% -NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES -% -CFL_NUMBER= 60 -CFL_REDUCTION_SPECIES= 1.0 -CFL_REDUCTION_TURB= 1.0 -% -ITER= 10 -% -% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% -% -LINEAR_SOLVER= FGMRES -LINEAR_SOLVER_PREC= ILU -LINEAR_SOLVER_ERROR= 1E-8 -LINEAR_SOLVER_ITER= 5 -% -% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% -% -CONV_NUM_METHOD_FLOW= FDS -MUSCL_FLOW= YES -SLOPE_LIMITER_FLOW = NONE -TIME_DISCRE_FLOW= EULER_IMPLICIT -% -% -------------------- SCALAR TRANSPORT ---------------------------------------% -% -KIND_SCALAR_MODEL= SPECIES_TRANSPORT -DIFFUSIVITY_MODEL= CONSTANT_DIFFUSIVITY -DIFFUSIVITY_CONSTANT= 0.001 -% -CONV_NUM_METHOD_SPECIES= BOUNDED_SCALAR -MUSCL_SPECIES= NO -SLOPE_LIMITER_SPECIES = NONE -% -TIME_DISCRE_SPECIES= EULER_IMPLICIT -% -SPECIES_INIT= 1.0 -SPECIES_CLIPPING= YES -SPECIES_CLIPPING_MIN= 0.0 -SPECIES_CLIPPING_MAX= 1.0 -% -% -------------------- TURBULENT TRANSPORT ---------------------------------------% -% -CONV_NUM_METHOD_TURB= BOUNDED_SCALAR -MUSCL_TURB= NO -% -% --------------------------- CONVERGENCE PARAMETERS --------------------------% -% -CONV_FIELD= RMS_PRESSURE, RMS_VELOCITY-X, RMS_VELOCITY-Y, RMS_TKE, RMS_SPECIES -CONV_RESIDUAL_MINVAL= -18 -CONV_STARTITER= 10 -% -% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% -% -MESH_FILENAME= primitiveVenturi.su2 -SCREEN_OUTPUT= INNER_ITER WALL_TIME \ - RMS_PRESSURE RMS_VELOCITY-X RMS_VELOCITY-Y RMS_TKE RMS_DISSIPATION RMS_SPECIES_0 \ - LINSOL_ITER LINSOL_RESIDUAL \ - LINSOL_ITER_TURB LINSOL_RESIDUAL_TURB \ - LINSOL_ITER_SPECIES LINSOL_RESIDUAL_SPECIES SURFACE_SPECIES_VARIANCE -SCREEN_WRT_FREQ_INNER= 10 -% -HISTORY_OUTPUT= RMS_RES FLOW_COEFF LINSOL SPECIES_COEFF SPECIES_COEFF_SURF -MARKER_ANALYZE= outlet gas_inlet air_axial_inlet -% -OUTPUT_FILES= RESTART_ASCII, PARAVIEW_MULTIBLOCK -VOLUME_OUTPUT= RESIDUAL, PRIMITIVE -OUTPUT_WRT_FREQ= 1000 -% -RESTART_SOL= NO -SOLUTION_FILENAME= solution -RESTART_FILENAME= restart -% -WRT_PERFORMANCE= YES diff --git a/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg b/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg index e914bf4b1aa2..b6104c898439 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg @@ -2,10 +2,10 @@ % % % SU2 configuration file % % Case description: Species mixing with 3 species, i.e. 2 transport equations % -% using inlet file for 2 inlets % -% Author: T. Kattmann % -% Institution: Bosch Thermotechniek B.V. % -% Date: 2021/10/14 % +% using FLUX and VALUE boundary conditions for the wall % +% Author: N. Beishuizen % +% Institution: TU Eindhoven % +% Date: 17-12-2025 % % File Version 8.3.0 "Harrier" % % % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -51,9 +51,10 @@ INLET_FILENAME= inlet_venturi.dat INC_INLET_TYPE= VELOCITY_INLET VELOCITY_INLET MARKER_INLET= ( gas_inlet, 300, 1.0, 1.0, 0.0, 0.0,\ air_axial_inlet, 300, 1.0, 0.0, -1.0, 0.0 ) -MARKER_INLET_SPECIES= (gas_inlet, 0.5, 0.5,\ - air_axial_inlet, 0.6, 0.0 ) -MARKER_WALL_SPECIES=(axis,FLUX,0.1,FLUX,0.2,wall,VALUE,0.5,VALUE,0.6) +MARKER_INLET_SPECIES= (gas_inlet, 0.0, 0.1,\ + air_axial_inlet, 0.1, 0.0 ) +% VALUE (Dirichlet) for species 1, FLUX (Neumann) for species 2 +MARKER_WALL_SPECIES=( wall, VALUE, 1.0, FLUX, 0.01 ) % INC_OUTLET_TYPE= PRESSURE_OUTLET MARKER_OUTLET= ( outlet, 0.0 ) @@ -90,7 +91,7 @@ KIND_SCALAR_MODEL= SPECIES_TRANSPORT DIFFUSIVITY_MODEL= CONSTANT_DIFFUSIVITY DIFFUSIVITY_CONSTANT= 0.001 % -CONV_NUM_METHOD_SPECIES= SCALAR_UPWIND +CONV_NUM_METHOD_SPECIES= BOUNDED_SCALAR MUSCL_SPECIES= NO SLOPE_LIMITER_SPECIES = NONE % @@ -103,13 +104,13 @@ SPECIES_CLIPPING_MAX= 1.0, 1.0 % % -------------------- TURBULENT TRANSPORT ---------------------------------------% % -CONV_NUM_METHOD_TURB= SCALAR_UPWIND +CONV_NUM_METHOD_TURB= BOUNDED_SCALAR MUSCL_TURB= NO % % --------------------------- CONVERGENCE PARAMETERS --------------------------% % CONV_FIELD= RMS_PRESSURE, RMS_VELOCITY-X, RMS_VELOCITY-Y, RMS_TKE, RMS_SPECIES -CONV_RESIDUAL_MINVAL= -18 +CONV_RESIDUAL_MINVAL= -16 CONV_STARTITER= 10 % % ------------------------- INPUT/OUTPUT INFORMATION --------------------------% diff --git a/config_template.cfg b/config_template.cfg index 97f5d967970b..9b008e2bccdd 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -1057,6 +1057,7 @@ MARKER_ISOTHERMAL= ( NONE ) % Format: ( marker_name, BC_TYPE_1, value_1, BC_TYPE_2, value_2, ..., BC_TYPE_N, value_N, ... ) % where BC_TYPE is either FLUX (Neumann) or VALUE (Dirichlet) for each species % Example for 3 species: MARKER_WALL_SPECIES= (wall1, FLUX, 0.0, VALUE, 0.3, VALUE, 0.2, wall2, VALUE, 0.1, VALUE, 0.4, FLUX, 0.0) +% If MARKER_WALL_SPECIES is not specified, a zero-flux boundary condition will be used. MARKER_WALL_SPECIES= ( NONE ) % % Far-field boundary marker(s) (NONE = no marker) From bbebbbb317f0829d9fd6ab5940b1dc7b65c65517 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Wed, 17 Dec 2025 11:08:20 +0100 Subject: [PATCH 14/31] empty line --- Common/src/CConfig.cpp | 1 - 1 file changed, 1 deletion(-) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index c0adf00428d8..06bba4bd025a 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -5739,7 +5739,6 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i /*--- Check whether some variables (or their sums) are in physical bounds. [0,1] for species related quantities. ---*/ /*--- Note, only for species transport, not for flamelet model ---*/ - if (Kind_Species_Model == SPECIES_MODEL::SPECIES_TRANSPORT) { su2double Species_Init_Sum = 0.0; for (unsigned short iSpecies = 0; iSpecies < nSpecies; iSpecies++) { From b916b035762b90a29dc8fc5997a25fa6c8325c3b Mon Sep 17 00:00:00 2001 From: bigfooted Date: Wed, 17 Dec 2025 11:20:45 +0100 Subject: [PATCH 15/31] add regression test --- TestCases/parallel_regression.py | 8 ++++++++ .../species3_primitiveVenturi_flux_value.cfg | 11 +++++------ 2 files changed, 13 insertions(+), 6 deletions(-) diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index c5a3a055ad5a..f0f903ce1e3e 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -1662,6 +1662,14 @@ def main(): species3_primitiveVenturi_inletFile.test_vals = [-5.537438, -4.503863, -4.553632, -5.400874, -0.945967, -5.818774, -5.945211, 5.000000, -0.544749, 5.000000, -2.599435, 5.000000, -0.596360] test_list.append(species3_primitiveVenturi_inletFile) + # 3 species (2 eq) primitive venturi mixing with new flux and value boundary conditions + species3_primitiveVenturi_fluxvalue = TestCase('species3_primitiveVenturi_fluxvalue') + species3_primitiveVenturi_fluxvalue.cfg_dir = "species_transport/venturi_primitive_3species" + species3_primitiveVenturi_fluxvalue.cfg_file = "species3_primitiveVenturi_flux_value.cfg" + species3_primitiveVenturi_fluxvalue.test_iter = 50 + species3_primitiveVenturi_fluxvalue.test_vals = [-5.537438, -4.503863, -4.553632, -5.400874, -0.945967, -5.818774, -5.945211, 5.000000, -0.544749, 5.000000, -2.599435, 5.000000, -0.596360] + test_list.append(species3_primitiveVenturi_inletFile) + # rectangle passive transport validation species_passive_val = TestCase('species_passive_val') species_passive_val.cfg_dir = "species_transport/passive_transport_validation" diff --git a/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg b/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg index b6104c898439..68880a6cbfca 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_flux_value.cfg @@ -75,7 +75,7 @@ ITER= 100 % LINEAR_SOLVER= FGMRES LINEAR_SOLVER_PREC= ILU -LINEAR_SOLVER_ERROR= 1E-8 +LINEAR_SOLVER_ERROR= 1E-3 LINEAR_SOLVER_ITER= 5 % % -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% @@ -118,9 +118,7 @@ CONV_STARTITER= 10 MESH_FILENAME= primitiveVenturi.su2 SCREEN_OUTPUT= INNER_ITER WALL_TIME \ RMS_PRESSURE RMS_VELOCITY-X RMS_VELOCITY-Y RMS_TKE RMS_DISSIPATION RMS_SPECIES_0 RMS_SPECIES_1 \ - LINSOL_ITER LINSOL_RESIDUAL \ - LINSOL_ITER_TURB LINSOL_RESIDUAL_TURB \ - LINSOL_ITER_SPECIES LINSOL_RESIDUAL_SPECIES + SURFACE_SPECIES_0 SURFACE_SPECIES_1 SCREEN_WRT_FREQ_INNER= 10 % HISTORY_OUTPUT= ITER RMS_RES LINSOL SPECIES_COEFF SPECIES_COEFF_SURF @@ -128,7 +126,8 @@ CONV_FILENAME= history_inlet MARKER_ANALYZE= gas_inlet, air_axial_inlet, outlet MARKER_ANALYZE_AVERAGE= AREA % -OUTPUT_FILES= RESTART_ASCII, PARAVIEW_MULTIBLOCK +%OUTPUT_FILES= RESTART_ASCII, PARAVIEW_MULTIBLOCK +OUTPUT_FILES= NONE VOLUME_OUTPUT= RESIDUAL, PRIMITIVE OUTPUT_WRT_FREQ= 1000 % @@ -136,4 +135,4 @@ RESTART_SOL= NO SOLUTION_FILENAME= solution RESTART_FILENAME= restart % -WRT_PERFORMANCE= YES +WRT_PERFORMANCE= NO From 776059ddea398e3da3fdf5f8d5e90de005ddb6a3 Mon Sep 17 00:00:00 2001 From: Nijso Date: Fri, 19 Dec 2025 23:24:41 +0100 Subject: [PATCH 16/31] Apply suggestion from @bigfooted --- TestCases/parallel_regression.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index f0f903ce1e3e..bc7b2c319270 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -1668,7 +1668,7 @@ def main(): species3_primitiveVenturi_fluxvalue.cfg_file = "species3_primitiveVenturi_flux_value.cfg" species3_primitiveVenturi_fluxvalue.test_iter = 50 species3_primitiveVenturi_fluxvalue.test_vals = [-5.537438, -4.503863, -4.553632, -5.400874, -0.945967, -5.818774, -5.945211, 5.000000, -0.544749, 5.000000, -2.599435, 5.000000, -0.596360] - test_list.append(species3_primitiveVenturi_inletFile) + test_list.append(species3_primitiveVenturi_fluxvalue) # rectangle passive transport validation species_passive_val = TestCase('species_passive_val') From d18433df5157019014454ae13b6e5bf9d1a19c66 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Tue, 23 Dec 2025 21:43:44 +0100 Subject: [PATCH 17/31] update regression --- TestCases/parallel_regression.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index bc7b2c319270..91b900836279 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -1667,7 +1667,7 @@ def main(): species3_primitiveVenturi_fluxvalue.cfg_dir = "species_transport/venturi_primitive_3species" species3_primitiveVenturi_fluxvalue.cfg_file = "species3_primitiveVenturi_flux_value.cfg" species3_primitiveVenturi_fluxvalue.test_iter = 50 - species3_primitiveVenturi_fluxvalue.test_vals = [-5.537438, -4.503863, -4.553632, -5.400874, -0.945967, -5.818774, -5.945211, 5.000000, -0.544749, 5.000000, -2.599435, 5.000000, -0.596360] + species3_primitiveVenturi_fluxvalue.test_vals = [-4.563229, -5.504499, -0.861681, -5.822963, -6.458352, 1.257908, 0.122218, 0.317705, 0.817985, 0.241494, 0.102507, 0.004981, 0.134006] test_list.append(species3_primitiveVenturi_fluxvalue) # rectangle passive transport validation From a9bc5dc2a5c323b6280dd546e8b280c07ba6f438 Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 24 Dec 2025 08:06:29 +0100 Subject: [PATCH 18/31] Update Common/include/option_structure.hpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- Common/include/option_structure.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Common/include/option_structure.hpp b/Common/include/option_structure.hpp index d8a560ab0fc4..8954b8c019a8 100644 --- a/Common/include/option_structure.hpp +++ b/Common/include/option_structure.hpp @@ -1838,7 +1838,7 @@ static const MapType Giles_Map = { /*! * \brief Types of wall species boundary conditions. */ -enum WALL_SPECIES_TYPE { +enum class WALL_SPECIES_TYPE { WALL_SPECIES_FLUX = 1, /*!< \brief Neumann flux boundary condition for wall species. */ WALL_SPECIES_VALUE = 2 /*!< \brief Dirichlet value boundary condition for wall species. */ }; From 51e84b508db1c73f172380d85716dc75897e8327 Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 24 Dec 2025 08:06:49 +0100 Subject: [PATCH 19/31] Update Common/include/CConfig.hpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- Common/include/CConfig.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 6819585b289e..de72b6702276 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -7162,7 +7162,7 @@ class CConfig { * \param[in] iSpecies - Species index. * \return The wall species type (WALL_SPECIES_FLUX or WALL_SPECIES_VALUE). */ - unsigned short GetWall_SpeciesType(const string& val_marker, unsigned short iSpecies) const; + WALL_SPECIES_TYPE GetWall_SpeciesType(const string& val_marker, unsigned short iSpecies) const; /*! * \brief Get the turbulent properties values at an inlet boundary From d46ac46bb422270407b70143a5d3a7751fa6f167 Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 24 Dec 2025 08:06:57 +0100 Subject: [PATCH 20/31] Update SU2_CFD/include/drivers/CDriverBase.hpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- SU2_CFD/include/drivers/CDriverBase.hpp | 1 - 1 file changed, 1 deletion(-) diff --git a/SU2_CFD/include/drivers/CDriverBase.hpp b/SU2_CFD/include/drivers/CDriverBase.hpp index 33b75fad0138..f05ab722372c 100644 --- a/SU2_CFD/include/drivers/CDriverBase.hpp +++ b/SU2_CFD/include/drivers/CDriverBase.hpp @@ -573,7 +573,6 @@ class CDriverBase { * \param[in] iVertex - Vertex identifier. * \param[in] WallScalar - Value of the normal heat flux. */ - inline void SetMarkerCustomScalar(unsigned short iMarker, unsigned long iVertex, vector WallScalar) { auto* solver = solver_container[selected_zone][INST_0][MESH_0][SPECIES_SOL]; solver->SetCustomBoundaryScalar(iMarker, iVertex, WallScalar); From 4aa52168dfca8543cacb0191bf990532c533e16e Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 24 Dec 2025 08:07:26 +0100 Subject: [PATCH 21/31] Update SU2_CFD/include/solvers/CSolver.hpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- SU2_CFD/include/solvers/CSolver.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/SU2_CFD/include/solvers/CSolver.hpp b/SU2_CFD/include/solvers/CSolver.hpp index 8f3f9c4136aa..fcbdb01bba20 100644 --- a/SU2_CFD/include/solvers/CSolver.hpp +++ b/SU2_CFD/include/solvers/CSolver.hpp @@ -2857,7 +2857,7 @@ class CSolver { unsigned short val_dim) const { return 0; } - /*! + /*! * \brief Set the value of the customized normal scalar values/flux at a specified vertex on a specified marker. * \param[in] val_marker - Marker value * \param[in] val_vertex - Boundary vertex value From 9494da4fd8cb6a1b978655ec6f33f41b2c46ffff Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 24 Dec 2025 08:07:43 +0100 Subject: [PATCH 22/31] Update SU2_CFD/include/solvers/CSpeciesSolver.hpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- SU2_CFD/include/solvers/CSpeciesSolver.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/SU2_CFD/include/solvers/CSpeciesSolver.hpp b/SU2_CFD/include/solvers/CSpeciesSolver.hpp index 8b15f6dfa9ad..cf82c121ad5d 100644 --- a/SU2_CFD/include/solvers/CSpeciesSolver.hpp +++ b/SU2_CFD/include/solvers/CSpeciesSolver.hpp @@ -241,7 +241,7 @@ class CSpeciesSolver : public CScalarSolver { * \param[in] val_customBoundaryScalar - Vector of scalar values */ inline void SetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, - vector val_customBoundaryScalar) { + vector val_customBoundaryScalar) final { for (auto iVar = 0u; iVar < nVar; iVar++) { CustomBoundaryScalar[val_marker](val_vertex, iVar) = val_customBoundaryScalar[iVar]; } From 76c604e919231973d89e853b38be50c62203d24f Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 24 Dec 2025 08:08:02 +0100 Subject: [PATCH 23/31] Update SU2_CFD/src/solvers/CSpeciesSolver.cpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- SU2_CFD/src/solvers/CSpeciesSolver.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/SU2_CFD/src/solvers/CSpeciesSolver.cpp b/SU2_CFD/src/solvers/CSpeciesSolver.cpp index 62154c735d4c..f26332ebe798 100644 --- a/SU2_CFD/src/solvers/CSpeciesSolver.cpp +++ b/SU2_CFD/src/solvers/CSpeciesSolver.cpp @@ -428,7 +428,7 @@ void CSpeciesSolver::BC_Isothermal_Wall(CGeometry* geometry, CSolver** solver_co void CSpeciesSolver::BC_HeatFlux_Wall(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, CNumerics* visc_numerics, CConfig* config, unsigned short val_marker) { - BC_Wall_Generic(geometry, solver_container, conv_numerics, visc_numerics, config, val_marker,0); + BC_Wall_Generic(geometry, solver_container, conv_numerics, visc_numerics, config, val_marker, false); } void CSpeciesSolver::BC_Wall_Generic(CGeometry* geometry, CSolver** solver_container, From f71f958a89d65088d8b00e655c5937441bd54e04 Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 24 Dec 2025 08:08:27 +0100 Subject: [PATCH 24/31] Update SU2_CFD/src/solvers/CSpeciesSolver.cpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- SU2_CFD/src/solvers/CSpeciesSolver.cpp | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/SU2_CFD/src/solvers/CSpeciesSolver.cpp b/SU2_CFD/src/solvers/CSpeciesSolver.cpp index f26332ebe798..d2e839a1a3a4 100644 --- a/SU2_CFD/src/solvers/CSpeciesSolver.cpp +++ b/SU2_CFD/src/solvers/CSpeciesSolver.cpp @@ -432,8 +432,7 @@ void CSpeciesSolver::BC_HeatFlux_Wall(CGeometry* geometry, CSolver** solver_cont } void CSpeciesSolver::BC_Wall_Generic(CGeometry* geometry, CSolver** solver_container, - CNumerics* conv_numerics, CNumerics* visc_numerics, - CConfig* config, unsigned short val_marker, bool cht_mode) { + CConfig* config, unsigned short val_marker) { const bool implicit = config->GetKind_TimeIntScheme_Flow() == EULER_IMPLICIT; const bool py_custom = config->GetMarker_All_PyCustom(val_marker); From ae628502c08efd7889464214f63736d2e2a0924d Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 24 Dec 2025 08:21:42 +0100 Subject: [PATCH 25/31] Update SU2_CFD/include/solvers/CSpeciesSolver.hpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- SU2_CFD/include/solvers/CSpeciesSolver.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/SU2_CFD/include/solvers/CSpeciesSolver.hpp b/SU2_CFD/include/solvers/CSpeciesSolver.hpp index cf82c121ad5d..89819337cc7a 100644 --- a/SU2_CFD/include/solvers/CSpeciesSolver.hpp +++ b/SU2_CFD/include/solvers/CSpeciesSolver.hpp @@ -41,7 +41,7 @@ class CSpeciesSolver : public CScalarSolver { unsigned short Inlet_Position; /*!< \brief Column index for scalar variables in inlet files. */ vector Inlet_SpeciesVars; /*!< \brief Species variables at inlet profiles. */ vector Wall_SpeciesVars; /*!< \brief Species variables at profiles. */ - vector >CustomBoundaryScalar; + vector> CustomBoundaryScalar; public: /*! From 9749e64ee96ee162e6fcf2501d3fde7539c2608c Mon Sep 17 00:00:00 2001 From: Nijso Date: Wed, 24 Dec 2025 08:25:11 +0100 Subject: [PATCH 26/31] Update SU2_CFD/include/solvers/CSpeciesSolver.hpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- SU2_CFD/include/solvers/CSpeciesSolver.hpp | 11 ----------- 1 file changed, 11 deletions(-) diff --git a/SU2_CFD/include/solvers/CSpeciesSolver.hpp b/SU2_CFD/include/solvers/CSpeciesSolver.hpp index 89819337cc7a..2feef978a5a1 100644 --- a/SU2_CFD/include/solvers/CSpeciesSolver.hpp +++ b/SU2_CFD/include/solvers/CSpeciesSolver.hpp @@ -247,15 +247,4 @@ class CSpeciesSolver : public CScalarSolver { } } - /*! - * \brief Get custom boundary scalar value for a specific variable. - * \param[in] val_marker - Boundary marker index - * \param[in] val_vertex - Boundary vertex index - * \param[in] iVar - Variable index - * \return Custom boundary scalar value - */ - inline su2double GetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, unsigned short iVar) const { - return CustomBoundaryScalar[val_marker](val_vertex, iVar); - } - }; From 72eb31e1248dc185cf88465dd5a4927aa2ace6e1 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Wed, 24 Dec 2025 08:25:47 +0100 Subject: [PATCH 27/31] implement review changes --- Common/src/CConfig.cpp | 4 +- .../include/solvers/CFVMFlowSolverBase.hpp | 8 --- SU2_CFD/include/solvers/CSpeciesSolver.hpp | 4 +- SU2_CFD/src/solvers/CSpeciesSolver.cpp | 65 +++++++++---------- 4 files changed, 34 insertions(+), 47 deletions(-) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 06bba4bd025a..7147bff8d976 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -9262,7 +9262,7 @@ su2double CConfig::GetWall_SpeciesVal(const string& val_marker, unsigned short i return 0.0; } -unsigned short CConfig::GetWall_SpeciesType(const string& val_marker, unsigned short iSpecies) const { +WALL_SPECIES_TYPE CConfig::GetWall_SpeciesType(const string& val_marker, unsigned short iSpecies) const { /*--- Search for the marker in the wall species list ---*/ for (unsigned short iMarker_Wall_Species = 0; iMarker_Wall_Species < nMarker_Wall_Species; iMarker_Wall_Species++) { if (Marker_Wall_Species[iMarker_Wall_Species] == val_marker) { @@ -9270,7 +9270,7 @@ unsigned short CConfig::GetWall_SpeciesType(const string& val_marker, unsigned s } } /*--- If marker not found (MARKER_WALL_SPECIES=NONE), return FLUX type (zero flux BC) ---*/ - return WALL_SPECIES_FLUX; + return WALL_SPECIES_TYPE::WALL_SPECIES_FLUX; } const su2double* CConfig::GetInlet_TurbVal(const string& val_marker) const { diff --git a/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp b/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp index 2fc3ae88f4f7..362e2d51c4b0 100644 --- a/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp +++ b/SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp @@ -2202,14 +2202,6 @@ class CFVMFlowSolverBase : public CSolver { Inlet_FlowDir[val_marker][val_vertex][val_dim] = val_flowdir; } - /*! - * \brief Set the value of the customized normal scalar values/flux at a specified vertex on a specified marker. - * \param[in] val_marker - Marker value - * \param[in] val_vertex - Boundary vertex value - */ - inline void SetCustomBoundaryScalar(unsigned short val_marker, unsigned long val_vertex, - vector val_customBoundaryScalar) { } - /*! * \brief Update the multi-grid structure for the customized boundary conditions. * \param geometry_container - Geometrical definition. diff --git a/SU2_CFD/include/solvers/CSpeciesSolver.hpp b/SU2_CFD/include/solvers/CSpeciesSolver.hpp index cf82c121ad5d..351a5d666d0b 100644 --- a/SU2_CFD/include/solvers/CSpeciesSolver.hpp +++ b/SU2_CFD/include/solvers/CSpeciesSolver.hpp @@ -149,7 +149,7 @@ class CSpeciesSolver : public CScalarSolver { CConfig* config, unsigned short val_marker) final; /*! - * \brief Impose the isothermal wall boundary condition. + * \brief Impose the isothermal wall Dirichlet boundary condition (value). * \param[in] geometry - Geometrical definition of the problem. * \param[in] solver_container - Container vector with all the solutions. * \param[in] conv_numerics - Description of the numerical method. @@ -161,7 +161,7 @@ class CSpeciesSolver : public CScalarSolver { CNumerics* visc_numerics, CConfig* config, unsigned short val_marker) override; /*! - * \brief Impose the heat flux wall boundary condition. + * \brief Impose the heat flux Neumann wall boundary condition (flux). * \param[in] geometry - Geometrical definition of the problem. * \param[in] solver_container - Container vector with all the solutions. * \param[in] conv_numerics - Description of the numerical method. diff --git a/SU2_CFD/src/solvers/CSpeciesSolver.cpp b/SU2_CFD/src/solvers/CSpeciesSolver.cpp index d2e839a1a3a4..677ef947ffa9 100644 --- a/SU2_CFD/src/solvers/CSpeciesSolver.cpp +++ b/SU2_CFD/src/solvers/CSpeciesSolver.cpp @@ -441,50 +441,45 @@ void CSpeciesSolver::BC_Wall_Generic(CGeometry* geometry, CSolver** solver_conta SU2_OMP_FOR_DYN(OMP_MIN_SIZE) for (auto iVertex = 0u; iVertex < geometry->nVertex[val_marker]; iVertex++) { - const auto iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); - - if (!geometry->nodes->GetDomain(iPoint)) continue; - - const auto Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); - - su2double Area = GeometryToolbox::Norm(nDim, Normal); + const auto iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); - for (auto iVar = 0u; iVar < nVar; iVar++) { + if (!geometry->nodes->GetDomain(iPoint)) continue; - // Get wall species boundary condition type and value for this marker and species - su2double WallSpeciesValue = config->GetWall_SpeciesVal(Marker_Tag, iVar); - unsigned short wallspeciestype = config->GetWall_SpeciesType(Marker_Tag, iVar); + const auto Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); - su2double WallSpecies = WallSpeciesValue; + su2double Area = GeometryToolbox::Norm(nDim, Normal); - /*--- Get the scalar values from the python wrapper. ---*/ - if (py_custom) { - WallSpecies = GetCustomBoundaryScalar(val_marker, iVertex, iVar); - } + for (auto iVar = 0u; iVar < nVar; iVar++) { - switch(wallspeciestype){ - case WALL_SPECIES_FLUX: - //Flux Boundary condition - LinSysRes(iPoint, iVar) -= WallSpecies * Area; + // Get wall species boundary condition type and value for this marker and species + su2double WallSpeciesValue = config->GetWall_SpeciesVal(Marker_Tag, iVar); + unsigned short wallspeciestype = config->GetWall_SpeciesType(Marker_Tag, iVar); - break; - case WALL_SPECIES_VALUE: - //Dirichlet Strong Boundary Condition - nodes->SetSolution(iPoint, iVar, WallSpecies); - nodes->SetSolution_Old(iPoint, iVar, WallSpecies); + su2double WallSpecies = WallSpeciesValue; - LinSysRes(iPoint, iVar) = 0.0; + /*--- Get the scalar values from the python wrapper. ---*/ + if (py_custom) { + WallSpecies = GetCustomBoundaryScalar(val_marker, iVertex, iVar); + } - if (implicit) { - unsigned long total_index = iPoint * nVar + iVar; - Jacobian.DeleteValsRowi(total_index); - } - break; + switch(wallspeciestype) { + case WALL_SPECIES_TYPE::WALL_SPECIES_FLUX: + //Flux Boundary condition + LinSysRes(iPoint, iVar) -= WallSpecies * Area; + break; + case WALL_SPECIES_TYPE::WALL_SPECIES_VALUE: + //Dirichlet Strong Boundary Condition + nodes->SetSolution(iPoint, iVar, WallSpecies); + nodes->SetSolution_Old(iPoint, iVar, WallSpecies); + LinSysRes(iPoint, iVar) = 0.0; + if (implicit) { + unsigned long total_index = iPoint * nVar + iVar; + Jacobian.DeleteValsRowi(total_index); + } + break; } - } - } - END_SU2_OMP_FOR - + } + END_SU2_OMP_FOR } From 3357529c0d1deedbf3deef346d32a6222365b803 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Wed, 24 Dec 2025 08:51:10 +0100 Subject: [PATCH 28/31] implement review changes --- Common/include/CConfig.hpp | 4 +- Common/include/option_structure.hpp | 8 +-- Common/include/option_structure.inl | 8 +-- Common/src/CConfig.cpp | 2 +- SU2_CFD/include/solvers/CSpeciesSolver.hpp | 17 +++--- SU2_CFD/src/solvers/CSpeciesSolver.cpp | 67 +++++++++++----------- 6 files changed, 53 insertions(+), 53 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index de72b6702276..1d49783cece9 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -297,7 +297,7 @@ class CConfig { su2double **Inlet_Velocity; /*!< \brief Specified flow velocity vectors for supersonic inlet boundaries. */ su2double **Inlet_SpeciesVal; /*!< \brief Specified species vector for inlet boundaries. */ su2double **Inlet_TurbVal; /*!< \brief Specified turbulent intensity and viscosity ratio for inlet boundaries. */ - unsigned short **Kind_Wall_Species; /*!< \brief Species boundary condition type for wall boundaries (FLUX or VALUE) per species. */ + WALL_SPECIES_TYPE **Kind_Wall_Species; /*!< \brief Species boundary condition type for wall boundaries (FLUX or VALUE) per species. */ su2double **Wall_SpeciesVal; /*!< \brief Specified species flux or value for wall boundaries per species. */ su2double *EngineInflow_Target; /*!< \brief Specified fan face targets for nacelle boundaries. */ su2double *Inflow_Mach; /*!< \brief Specified fan face mach for nacelle boundaries. */ @@ -1399,7 +1399,7 @@ class CConfig { template void addWallSpeciesOption(const string name, unsigned short & nMarker_Wall_Species, string * & Marker_Wall_Species, - unsigned short** & option_field, const map & enum_map, + WALL_SPECIES_TYPE** & option_field, const map & enum_map, su2double** & value, unsigned short & nSpecies_per_Wall); void addExhaustOption(const string& name, unsigned short & nMarker_Exhaust, string * & Marker_Exhaust, diff --git a/Common/include/option_structure.hpp b/Common/include/option_structure.hpp index 8954b8c019a8..2ec2454fba09 100644 --- a/Common/include/option_structure.hpp +++ b/Common/include/option_structure.hpp @@ -1839,12 +1839,12 @@ static const MapType Giles_Map = { * \brief Types of wall species boundary conditions. */ enum class WALL_SPECIES_TYPE { - WALL_SPECIES_FLUX = 1, /*!< \brief Neumann flux boundary condition for wall species. */ - WALL_SPECIES_VALUE = 2 /*!< \brief Dirichlet value boundary condition for wall species. */ + WALL_SPECIES_FLUX, /*!< \brief Neumann flux boundary condition for wall species. */ + WALL_SPECIES_VALUE /*!< \brief Dirichlet value boundary condition for wall species. */ }; static const MapType Wall_Map = { - MakePair("FLUX", WALL_SPECIES_FLUX) - MakePair("VALUE", WALL_SPECIES_VALUE) + MakePair("FLUX", WALL_SPECIES_TYPE::WALL_SPECIES_FLUX) + MakePair("VALUE", WALL_SPECIES_TYPE::WALL_SPECIES_VALUE) }; /*! diff --git a/Common/include/option_structure.inl b/Common/include/option_structure.inl index 2f149d3b97c6..022993dda169 100644 --- a/Common/include/option_structure.inl +++ b/Common/include/option_structure.inl @@ -1309,13 +1309,13 @@ class COptionWallSpecies : public COptionBase { string name; // identifier for the option unsigned short& size; string*& marker; - unsigned short**& field; // Reference to the field name (now 2D: marker x species) + WALL_SPECIES_TYPE**& field; // Reference to the field name (now 2D: marker x species) su2double**& value; // Now 2D: marker x species unsigned short& nSpecies_per_Wall; public: COptionWallSpecies(string option_field_name, unsigned short& nMarker_Wall_Species, string*& Marker_Wall_Species, - unsigned short**& option_field, const map m, su2double**& value, + WALL_SPECIES_TYPE**& option_field, const map m, su2double**& value, unsigned short& nSpecies_per_Wall) : size(nMarker_Wall_Species), marker(Marker_Wall_Species), @@ -1416,11 +1416,11 @@ class COptionWallSpecies : public COptionBase { // Allocate arrays this->size = marker_indices.size(); this->marker = new string[this->size]; - this->field = new unsigned short*[this->size]; + this->field = new WALL_SPECIES_TYPE*[this->size]; this->value = new su2double*[this->size]; for (unsigned short i = 0; i < this->size; i++) { - this->field[i] = new unsigned short[this->nSpecies_per_Wall]; + this->field[i] = new WALL_SPECIES_TYPE[this->nSpecies_per_Wall]; this->value[i] = new su2double[this->nSpecies_per_Wall]; } diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 7147bff8d976..2ac9c26568ba 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -515,7 +515,7 @@ void CConfig::addRiemannOption(const string name, unsigned short & nMarker_Riema template void CConfig::addWallSpeciesOption(const string name, unsigned short & nMarker_Wall_Species, string * & Marker_Wall_Species, - unsigned short** & option_field, const map & enum_map, + WALL_SPECIES_TYPE** & option_field, const map & enum_map, su2double** & value, unsigned short & nSpecies_per_Wall) { assert(option_map.find(name) == option_map.end()); all_options.insert(pair(name, true)); diff --git a/SU2_CFD/include/solvers/CSpeciesSolver.hpp b/SU2_CFD/include/solvers/CSpeciesSolver.hpp index 66215d866fc5..011f1daddb5b 100644 --- a/SU2_CFD/include/solvers/CSpeciesSolver.hpp +++ b/SU2_CFD/include/solvers/CSpeciesSolver.hpp @@ -157,8 +157,9 @@ class CSpeciesSolver : public CScalarSolver { * \param[in] config - Definition of the particular problem. * \param[in] val_marker - Surface marker where the boundary condition is applied. */ - void BC_Isothermal_Wall(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, - CNumerics* visc_numerics, CConfig* config, unsigned short val_marker) override; + void BC_Isothermal_Wall(CGeometry* geometry, CSolver** solver_container, + CNumerics* conv_numerics, CNumerics* visc_numerics, + CConfig* config, unsigned short val_marker) override; /*! * \brief Impose the heat flux Neumann wall boundary condition (flux). @@ -169,21 +170,19 @@ class CSpeciesSolver : public CScalarSolver { * \param[in] config - Definition of the particular problem. * \param[in] val_marker - Surface marker where the boundary condition is applied. */ - void BC_HeatFlux_Wall(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, - CNumerics* visc_numerics, CConfig* config, unsigned short val_marker) override; + void BC_HeatFlux_Wall(CGeometry* geometry, CSolver** solver_container, + CNumerics* conv_numerics, CNumerics* visc_numerics, + CConfig* config, unsigned short val_marker) override; /*! * \brief Generic wall boundary condition implementation. * \param[in] geometry - Geometrical definition of the problem. * \param[in] solver_container - Container vector with all the solutions. - * \param[in] conv_numerics - Description of the numerical method. - * \param[in] visc_numerics - Description of the numerical method. * \param[in] config - Definition of the particular problem. * \param[in] val_marker - Surface marker where the boundary condition is applied. - * \param[in] cht_mode - Flag for conjugate heat transfer mode (default: false). */ - void BC_Wall_Generic(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, - CNumerics* visc_numerics, CConfig* config, unsigned short val_marker, bool cht_mode = false); + void BC_Wall_Generic(CGeometry* geometry, CSolver** solver_container, + CConfig* config, unsigned short val_marker); /*--- Note that BC_Sym_Plane, BC_HeatFlux_Wall, BC_Isothermal_Wall are all treated as zero-flux BC for the * mass-factions, which can be fulfilled by no additional residual contribution on these nodes. diff --git a/SU2_CFD/src/solvers/CSpeciesSolver.cpp b/SU2_CFD/src/solvers/CSpeciesSolver.cpp index 677ef947ffa9..e4d21fcb1d69 100644 --- a/SU2_CFD/src/solvers/CSpeciesSolver.cpp +++ b/SU2_CFD/src/solvers/CSpeciesSolver.cpp @@ -420,15 +420,15 @@ void CSpeciesSolver::BC_Inlet(CGeometry* geometry, CSolver** solver_container, C void CSpeciesSolver::BC_Isothermal_Wall(CGeometry* geometry, CSolver** solver_container, - CNumerics* conv_numerics, CNumerics* visc_numerics, CConfig* config, - unsigned short val_marker) { - BC_Wall_Generic(geometry, solver_container, conv_numerics, visc_numerics, config, val_marker); + CNumerics* conv_numerics, CNumerics* visc_numerics, + CConfig* config, unsigned short val_marker) { + BC_Wall_Generic(geometry, solver_container, config, val_marker); } void CSpeciesSolver::BC_HeatFlux_Wall(CGeometry* geometry, CSolver** solver_container, - CNumerics* conv_numerics, CNumerics* visc_numerics, CConfig* config, - unsigned short val_marker) { - BC_Wall_Generic(geometry, solver_container, conv_numerics, visc_numerics, config, val_marker, false); + CNumerics* conv_numerics, CNumerics* visc_numerics, + CConfig* config, unsigned short val_marker) { + BC_Wall_Generic(geometry, solver_container, config, val_marker); } void CSpeciesSolver::BC_Wall_Generic(CGeometry* geometry, CSolver** solver_container, @@ -438,48 +438,49 @@ void CSpeciesSolver::BC_Wall_Generic(CGeometry* geometry, CSolver** solver_conta string Marker_Tag = config->GetMarker_All_TagBound(val_marker); - SU2_OMP_FOR_DYN(OMP_MIN_SIZE) - for (auto iVertex = 0u; iVertex < geometry->nVertex[val_marker]; iVertex++) { + for (auto iVar = 0u; iVar < nVar; iVar++) { - const auto iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); + // Get wall species boundary condition type and value for this marker and species + const su2double WallSpeciesValue = config->GetWall_SpeciesVal(Marker_Tag, iVar); + const WALL_SPECIES_TYPE wallspeciestype = config->GetWall_SpeciesType(Marker_Tag, iVar); - if (!geometry->nodes->GetDomain(iPoint)) continue; + SU2_OMP_FOR_DYN(OMP_MIN_SIZE) + for (auto iVertex = 0u; iVertex < geometry->nVertex[val_marker]; iVertex++) { - const auto Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); + const auto iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); - su2double Area = GeometryToolbox::Norm(nDim, Normal); + if (!geometry->nodes->GetDomain(iPoint)) continue; - for (auto iVar = 0u; iVar < nVar; iVar++) { + const auto Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); - // Get wall species boundary condition type and value for this marker and species - su2double WallSpeciesValue = config->GetWall_SpeciesVal(Marker_Tag, iVar); - unsigned short wallspeciestype = config->GetWall_SpeciesType(Marker_Tag, iVar); + su2double Area = GeometryToolbox::Norm(nDim, Normal); su2double WallSpecies = WallSpeciesValue; /*--- Get the scalar values from the python wrapper. ---*/ if (py_custom) { - WallSpecies = GetCustomBoundaryScalar(val_marker, iVertex, iVar); + WallSpecies = CustomBoundaryScalar[val_marker](iVertex,iVar); } - switch(wallspeciestype) { - case WALL_SPECIES_TYPE::WALL_SPECIES_FLUX: - //Flux Boundary condition - LinSysRes(iPoint, iVar) -= WallSpecies * Area; - break; - case WALL_SPECIES_TYPE::WALL_SPECIES_VALUE: - //Dirichlet Strong Boundary Condition - nodes->SetSolution(iPoint, iVar, WallSpecies); - nodes->SetSolution_Old(iPoint, iVar, WallSpecies); - LinSysRes(iPoint, iVar) = 0.0; - if (implicit) { - unsigned long total_index = iPoint * nVar + iVar; - Jacobian.DeleteValsRowi(total_index); - } - break; + switch(wallspeciestype) { + case WALL_SPECIES_TYPE::WALL_SPECIES_FLUX: + //Flux Boundary condition + LinSysRes(iPoint, iVar) -= WallSpecies * Area; + break; + case WALL_SPECIES_TYPE::WALL_SPECIES_VALUE: + //Dirichlet Strong Boundary Condition + nodes->SetSolution(iPoint, iVar, WallSpecies); + nodes->SetSolution_Old(iPoint, iVar, WallSpecies); + LinSysRes(iPoint, iVar) = 0.0; + if (implicit) { + unsigned long total_index = iPoint * nVar + iVar; + Jacobian.DeleteValsRowi(total_index); + } + break; + } } + END_SU2_OMP_FOR } - END_SU2_OMP_FOR } From f9b53b6c62bdb484475d00d40f46427328e296a1 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Thu, 25 Dec 2025 23:06:15 +0100 Subject: [PATCH 29/31] precommit --- Common/include/option_structure.inl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Common/include/option_structure.inl b/Common/include/option_structure.inl index 022993dda169..9dc2c829509f 100644 --- a/Common/include/option_structure.inl +++ b/Common/include/option_structure.inl @@ -1310,7 +1310,7 @@ class COptionWallSpecies : public COptionBase { unsigned short& size; string*& marker; WALL_SPECIES_TYPE**& field; // Reference to the field name (now 2D: marker x species) - su2double**& value; // Now 2D: marker x species + su2double**& value; // Now 2D: marker x species unsigned short& nSpecies_per_Wall; public: From 8828523759a2badc54c1b9d06c71da3dc2fd08bc Mon Sep 17 00:00:00 2001 From: Nijso Date: Mon, 29 Dec 2025 09:45:22 +0100 Subject: [PATCH 30/31] Update Common/include/option_structure.hpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- Common/include/option_structure.hpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Common/include/option_structure.hpp b/Common/include/option_structure.hpp index 23ac994d0683..0f5604ca69d5 100644 --- a/Common/include/option_structure.hpp +++ b/Common/include/option_structure.hpp @@ -1839,8 +1839,8 @@ static const MapType Giles_Map = { * \brief Types of wall species boundary conditions. */ enum class WALL_SPECIES_TYPE { - WALL_SPECIES_FLUX, /*!< \brief Neumann flux boundary condition for wall species. */ - WALL_SPECIES_VALUE /*!< \brief Dirichlet value boundary condition for wall species. */ + FLUX, /*!< \brief Neumann flux boundary condition for wall species. */ + VALUE /*!< \brief Dirichlet value boundary condition for wall species. */ }; static const MapType Wall_Map = { MakePair("FLUX", WALL_SPECIES_TYPE::WALL_SPECIES_FLUX) From 275afe3379ae4c05568ab610247eb4e908b00551 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Mon, 29 Dec 2025 11:20:06 +0100 Subject: [PATCH 31/31] change names --- Common/include/CConfig.hpp | 2 +- Common/include/option_structure.hpp | 4 ++-- Common/src/CConfig.cpp | 2 +- SU2_CFD/src/solvers/CSpeciesSolver.cpp | 4 ++-- 4 files changed, 6 insertions(+), 6 deletions(-) diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 8ac6961e655a..726495c42250 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -7171,7 +7171,7 @@ class CConfig { * \brief Get the species boundary condition type at a wall boundary for a specific species * \param[in] val_marker - Marker tag corresponding to the wall boundary. * \param[in] iSpecies - Species index. - * \return The wall species type (WALL_SPECIES_FLUX or WALL_SPECIES_VALUE). + * \return The wall species type (FLUX or VALUE). */ WALL_SPECIES_TYPE GetWall_SpeciesType(const string& val_marker, unsigned short iSpecies) const; diff --git a/Common/include/option_structure.hpp b/Common/include/option_structure.hpp index e9dcf43a4a3a..8375668778c4 100644 --- a/Common/include/option_structure.hpp +++ b/Common/include/option_structure.hpp @@ -1843,8 +1843,8 @@ enum class WALL_SPECIES_TYPE { VALUE /*!< \brief Dirichlet value boundary condition for wall species. */ }; static const MapType Wall_Map = { - MakePair("FLUX", WALL_SPECIES_TYPE::WALL_SPECIES_FLUX) - MakePair("VALUE", WALL_SPECIES_TYPE::WALL_SPECIES_VALUE) + MakePair("FLUX", WALL_SPECIES_TYPE::FLUX) + MakePair("VALUE", WALL_SPECIES_TYPE::VALUE) }; /*! diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index d0b086d51555..a608f288f7b3 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -9286,7 +9286,7 @@ WALL_SPECIES_TYPE CConfig::GetWall_SpeciesType(const string& val_marker, unsigne } } /*--- If marker not found (MARKER_WALL_SPECIES=NONE), return FLUX type (zero flux BC) ---*/ - return WALL_SPECIES_TYPE::WALL_SPECIES_FLUX; + return WALL_SPECIES_TYPE::FLUX; } const su2double* CConfig::GetInlet_TurbVal(const string& val_marker) const { diff --git a/SU2_CFD/src/solvers/CSpeciesSolver.cpp b/SU2_CFD/src/solvers/CSpeciesSolver.cpp index e4d21fcb1d69..309550740749 100644 --- a/SU2_CFD/src/solvers/CSpeciesSolver.cpp +++ b/SU2_CFD/src/solvers/CSpeciesSolver.cpp @@ -463,11 +463,11 @@ void CSpeciesSolver::BC_Wall_Generic(CGeometry* geometry, CSolver** solver_conta } switch(wallspeciestype) { - case WALL_SPECIES_TYPE::WALL_SPECIES_FLUX: + case WALL_SPECIES_TYPE::FLUX: //Flux Boundary condition LinSysRes(iPoint, iVar) -= WallSpecies * Area; break; - case WALL_SPECIES_TYPE::WALL_SPECIES_VALUE: + case WALL_SPECIES_TYPE::VALUE: //Dirichlet Strong Boundary Condition nodes->SetSolution(iPoint, iVar, WallSpecies); nodes->SetSolution_Old(iPoint, iVar, WallSpecies);