diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index d1ac360a4bb..54e00341482 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -697,6 +697,7 @@ class CConfig { nMarker_PyCustom, /*!< \brief Number of markers that are customizable in Python. */ nMarker_DV, /*!< \brief Number of markers affected by the design variables. */ nMarker_WallFunctions, /*!< \brief Number of markers for which wall functions must be applied. */ + nMarker_StrongBC, /*!< \brief Number of markers for which a strong BC must be applied. */ nMarker_SobolevBC; /*!< \brief Number of markers treaded in the gradient problem. */ string *Marker_Monitoring, /*!< \brief Markers to monitor. */ *Marker_Designing, /*!< \brief Markers to design. */ @@ -708,6 +709,7 @@ class CConfig { *Marker_PyCustom, /*!< \brief Markers that are customizable in Python. */ *Marker_DV, /*!< \brief Markers affected by the design variables. */ *Marker_WallFunctions, /*!< \brief Markers for which wall functions must be applied. */ + *Marker_StrongBC, /*!< \brief Markers for which a strong BC must be applied. */ *Marker_SobolevBC; /*!< \brief Markers in the gradient solver */ unsigned short nConfig_Files; /*!< \brief Number of config files for multiphysics problems. */ @@ -1179,7 +1181,7 @@ class CConfig { unsigned short maxBasisDim, /*!< \brief Maximum number of POD basis dimensions. */ rom_save_freq; /*!< \brief Frequency of unsteady time steps to save. */ - unsigned short nSpecies; /*!< \brief Number of transported species equations (for NEMO and species transport)*/ + unsigned short nSpecies = 0; /*!< \brief Number of transported species equations (for NEMO and species transport)*/ /* other NEMO configure options*/ unsigned short nSpecies_Cat_Wall, /*!< \brief No. of species for a catalytic wall. */ @@ -1210,6 +1212,23 @@ class CConfig { su2double* Species_Init; /*!< \brief Initial uniform value for scalar transport. */ unsigned short nSpecies_Init; /*!< \brief Number of entries of SPECIES_INIT */ + /*--- Additional flamelet solver options ---*/ + su2double flame_init[8]; /*!< \brief Initial solution parameters for flamelet solver.*/ + + /*--- lookup table ---*/ + unsigned short n_scalars = 0; /*!< \brief Number of transported scalars for flamelet LUT approach. */ + unsigned short n_lookups = 0; /*!< \brief Number of lookup variables, for visualization only. */ + unsigned short n_table_sources = 0; /*!< \brief Number of transported scalar source terms for LUT. */ + unsigned short n_user_scalars = 0; /*!< \brief Number of user defined (auxiliary) scalar transport equations. */ + unsigned short n_user_sources = 0; /*!< \brief Number of source terms for user defined (auxiliary) scalar transport equations. */ + unsigned short n_control_vars = 0; /*!< \brief Number of controlling variables (independent variables) for the LUT. */ + + vector table_scalar_names; /*!< \brief Names of transported scalar variables. */ + string* table_lookup_names; /*!< \brief Names of LUT variables. */ + string file_name_lut; /*!< \brief Filename of the LUT. */ + string* user_scalar_names; /*!< \brief Names of the user defined (auxiliary) transported scalars .*/ + string* user_source_names; /*!< \brief Names of the source terms for the user defined transported scalars. */ + /*! * \brief Set the default values of config options not set in the config file using another config object. * \param config - Config object to use the default values from. @@ -2101,6 +2120,63 @@ class CConfig { */ bool GetSpecies_StrongBC() const { return Species_StrongBC; } + /*! + * \brief Get the flame initialization. + * (x1,x2,x3) = flame offset. + * (x4,x5,x6) = flame normal, separating unburnt from burnt. + * (x7) = flame thickness, the length from unburnt to burnt conditions. + * (x8) = flame burnt thickness, the length to stay at burnt conditions. + * \return Flame initialization for the flamelet model. + */ + const su2double* GetFlameInit() const { return flame_init; } + + /*! + * \brief Get the number of control variables for flamelet model. + */ + unsigned short GetNControlVars() const { return n_control_vars; } + + /*! + * \brief Get the number of total transported scalars for flamelet model. + */ + unsigned short GetNScalars() const { return n_scalars; } + + /*! + * \brief Get the number of user scalars for flamelet model. + */ + unsigned short GetNUserScalars() const { return n_user_scalars; } + + /*! + * \brief Get the name of the user scalar. + */ + const string& GetUserScalarName(unsigned short i_user_scalar) const { + static const std::string none = "NONE"; + if (n_user_scalars > 0) return user_scalar_names[i_user_scalar]; else return none; + } + + /*! + * \brief Get the name of the user scalar source term. + */ + const string& GetUserSourceName(unsigned short i_user_source) const { + static const std::string none = "NONE"; + if (n_user_sources > 0) return user_source_names[i_user_source]; else return none; + } + + /*! + * \brief Get the number of transported scalars for combustion. + */ + unsigned short GetNLookups() const { return n_lookups; } + + /*! + * \brief Get the name of the variable that we want to retrieve from the lookup table. + */ + const string& GetLUTLookupName(unsigned short i_lookup) const { return table_lookup_names[i_lookup]; } + + /*! + * \brief Get the file name of the look up table. + * \return File name of the look up table. + */ + const string& GetFileNameLUT() const { return file_name_lut; } + /*! * \brief Get the Young's modulus of elasticity. * \return Value of the Young's modulus of elasticity. @@ -9305,6 +9381,13 @@ class CConfig { */ su2double GetWall_Emissivity(const string& val_index) const; + /*! + * \brief Get if boundary is strong or weak. + * \param[in] val_index - Index corresponding to the boundary. + * \return true if strong BC. + */ + bool GetMarker_StrongBC(const string& val_index) const; + /*! * \brief Get the value of the CFL condition for radiation solvers. * \return Value of the CFL condition for radiation solvers. diff --git a/Common/include/containers/CLookUpTable.hpp b/Common/include/containers/CLookUpTable.hpp index f6cd4c91bb8..0effac26d4d 100644 --- a/Common/include/containers/CLookUpTable.hpp +++ b/Common/include/containers/CLookUpTable.hpp @@ -60,13 +60,13 @@ class CLookUpTable { unsigned short table_dim = 2; /*!< \brief Table dimension.*/ /*! - * \brief the lower and upper limits of the z, y and x variable for each table level. + * \brief The lower and upper limits of the z, y and x variable for each table level. */ std::pair limits_table_z; su2vector> limits_table_y, limits_table_x; /*! \brief Holds the variable names stored in the table file. - * Order is in sync with data + * Order is in sync with data. */ su2vector names_var; @@ -99,12 +99,12 @@ class CLookUpTable { su2vector trap_map_x_y; /*! \brief - * vector of all the weight factors for the interpolation. + * Vector of all the weight factors for the interpolation. */ su2vector> interp_mat_inv_x_y; /*! \brief - * returns the index to the variable in the lookup table. + * Returns the index to the variable in the lookup table. */ inline unsigned int GetIndexOfVar(const std::string& nameVar) const { int index = find(names_var.begin(), names_var.end(), nameVar) - names_var.begin(); @@ -118,72 +118,84 @@ class CLookUpTable { return index; } + /*! \brief + * Returns true if the string is null or zero (ignores case). + */ + inline bool noSource(const std::string& name_var) const { + if (name_var.compare("NULL") == 0 || name_var.compare("Null") == 0 || name_var.compare("null") == 0 || + name_var.compare("ZERO") == 0 || name_var.compare("Zero") == 0 || name_var.compare("zero") == 0) { + return true; + } else { + return false; + } + } + /*! * \brief Get the pointer to the column data of the table (density, temperature, source terms, ...). - * \returns pointer to the column data. + * \returns Pointer to the column data. */ inline const su2double* GetDataP(const std::string& name_var, unsigned long i_level = 0) const { return table_data[i_level][GetIndexOfVar(name_var)]; } /*! - * \brief find the table limits, i.e. the minimum and maximum values of the 2 independent. + * \brief Find the table limits, i.e. the minimum and maximum values of the 2 independent * controlling variables. We put the values in the variables. * limits_table_x[2] and limit_table_y[2]. - * \param[in] name_CV1 - the string name for the first controlling variable. - * \param[in] name_CV2 - the string name of the second controlling variable. + * \param[in] name_CV1 - The string name for the first controlling variable. + * \param[in] name_CV2 - The string name of the second controlling variable. */ void FindTableLimits(const std::string& name_CV1, const std::string& name_CV2); /*! - * \brief construct a list of all the edges and a list of the pair of elements left and right of the edge. + * \brief Construct a list of all the edges and a list of the pair of elements left and right of the edge. */ void IdentifyUniqueEdges(); /*! - * \brief read the lookup table from file and store the data. + * \brief Read the lookup table from file and store the data. * \param[in] file_name_lut - the filename of the lookup table. */ void LoadTableRaw(const std::string& file_name_lut); /*! - * \brief compute vector of all (inverse) interpolation coefficients "interp_mat_inv_x_y" of all triangles. + * \brief Compute vector of all (inverse) interpolation coefficients "interp_mat_inv_x_y" of all triangles. */ void ComputeInterpCoeffs(); /*! - * \brief compute the inverse matrix for interpolation. - * \param[in] vec_CV1 - pointer to first coordinate (progress variable). - * \param[in] vec_CV2 - pointer to second coordinate (enthalpy). - * \param[in] point_ids - single triangle data. - * \param[out] interp_mat_inv - inverse matrix for interpolation. + * \brief Compute the inverse matrix for interpolation. + * \param[in] vec_CV1 - Pointer to first coordinate (progress variable). + * \param[in] vec_CV2 - Pointer to second coordinate (enthalpy). + * \param[in] point_ids - Single triangle data. + * \param[out] interp_mat_inv - Inverse matrix for interpolation. */ void GetInterpMatInv(const su2double* vec_CV1, const su2double* vec_CV2, std::array& point_ids, su2activematrix& interp_mat_inv); /*! - * \brief compute the interpolation coefficients for the triangular interpolation. - * \param[in] val_CV1 - value of first coordinate (progress variable). - * \param[in] val_CV2 - value of second coordinate (enthalpy). - * \param[in] interp_mat_inv - inverse matrix for interpolation. - * \param[out] interp_coeffs - interpolation coefficients. + * \brief Compute the interpolation coefficients for the triangular interpolation. + * \param[in] val_CV1 - Value of first coordinate (progress variable). + * \param[in] val_CV2 - Value of second coordinate (enthalpy). + * \param[in] interp_mat_inv - Inverse matrix for interpolation. + * \param[out] interp_coeffs - Interpolation coefficients. */ void GetInterpCoeffs(su2double val_CV1, su2double val_CV2, su2activematrix& interp_mat_inv, std::array& interp_coeffs); /*! - * \brief compute interpolated value of a point P in the triangle. - * \param[in] val_samples - pointer to the variable data. + * \brief Compute interpolated value of a point P in the triangle. + * \param[in] val_samples - Pointer to the variable data. * \param[in] val_triangle - ID to the triangle. - * \param[in] val_interp_coeffs - interpolation coefficients using the point data in P. - * \returns resulting value of the interpolation. + * \param[in] val_interp_coeffs - Interpolation coefficients using the point data in P. + * \returns Resulting value of the interpolation. */ su2double Interpolate(const su2double* val_samples, std::array& val_triangle, std::array& val_interp_coeffs); /*! * \brief Perform linear interpolation between two table levels for a single variable. - * \param[in] val_CV3 - value of the third controlling variable at the query point. + * \param[in] val_CV3 - Value of the third controlling variable at the query point. * \param[in] lower_level - Table level index of the table level directly below the query point. * \param[in] upper_level - Table level index of the table level directly above the query point. * \param[in] lower_value - Result from x-y interpolation on the lower table level. @@ -195,7 +207,7 @@ class CLookUpTable { /*! * \brief Perform linear interpolation between two table levels for a vector of variables. - * \param[in] val_CV3 - value of the third controlling variable at the query point. + * \param[in] val_CV3 - Value of the third controlling variable at the query point. * \param[in] lower_level - Table level index of the table level directly below the query point. * \param[in] upper_level - Table level index of the table level directly above the query point. * \param[in] lower_values - Results from x-y interpolation on the lower table level. @@ -207,52 +219,52 @@ class CLookUpTable { std::vector& var_vals) const; /*! - * \brief find the point on the hull (boundary of the table) that is closest to the point P(val_CV1,val_CV2). - * \param[in] val_CV1 - first coordinate of point P(val_CV1,val_CV2) to check. - * \param[in] val_CV2 - second coordinate of point P(val_CV1,val_CV2) to check. - * \returns point id of the nearest neighbor on the hull. + * \brief Find the point on the hull (boundary of the table) that is closest to the point P(val_CV1,val_CV2). + * \param[in] val_CV1 - First coordinate of point P(val_CV1,val_CV2) to check. + * \param[in] val_CV2 - Second coordinate of point P(val_CV1,val_CV2) to check. + * \returns Point id of the nearest neighbor on the hull. */ unsigned long FindNearestNeighborOnHull(su2double val_CV1, su2double val_CV2, unsigned long i_level = 0); /*! * \brief Interpolate data based on distance-weighted averaging on the nearest two table nodes. - * \param[in] val_CV1 - - first coordinate of point P(val_CV1,val_CV2) to check. - * \param[in] val_CV2 - second coordinate of point P(val_CV1,val_CV2) to check. - * \param[in] val_names_var - vector of string names of the variables to look up. - * \param[out] val_vars - pointer to the vector of stored values of the variables to look up. + * \param[in] val_CV1 - First coordinate of point P(val_CV1,val_CV2) to check. + * \param[in] val_CV2 - Second coordinate of point P(val_CV1,val_CV2) to check. + * \param[in] val_names_var - Vector of string names of the variables to look up. + * \param[out] val_vars - Pointer to the vector of stored values of the variables to look up. */ void InterpolateToNearestNeighbors(const su2double val_CV1, const su2double val_CV2, const std::vector& names_var, std::vector& var_vals, const unsigned long i_level = 0); /*! - * \brief Interpolate data based on distance-weighted averaging on the nearest two table nodes for a single variable - * \param[in] val_CV1 - - first coordinate of point P(val_CV1,val_CV2) to check. - * \param[in] val_CV2 - second coordinate of point P(val_CV1,val_CV2) to check. - * \param[in] name_var - variable name to look up. - * \param[out] var_val - pointer to the to be interpolated variable. + * \brief Interpolate data based on distance-weighted averaging on the nearest two table nodes for a single variable. + * \param[in] val_CV1 - First coordinate of point P(val_CV1,val_CV2) to check. + * \param[in] val_CV2 - Second coordinate of point P(val_CV1,val_CV2) to check. + * \param[in] name_var - Variable name to look up. + * \param[out] var_val - Pointer to the to be interpolated variable. */ void InterpolateToNearestNeighbors(const su2double val_CV1, const su2double val_CV2, const std::string& name_var, su2double* var_val, const unsigned long i_level = 0); /*! - * \brief determine if a point P(val_CV1,val_CV2) is inside the triangle val_id_triangle. - * \param[in] val_CV1 - first coordinate of point P(val_CV1,val_CV2) to check. - * \param[in] val_CV2 - second coordinate of point P(val_CV1,val_CV2) to check. + * \brief Determine if a point P(val_CV1,val_CV2) is inside the triangle val_id_triangle. + * \param[in] val_CV1 - First coordinate of point P(val_CV1,val_CV2) to check. + * \param[in] val_CV2 - Second coordinate of point P(val_CV1,val_CV2) to check. * \param[in] val_id_triangle - ID of the triangle to check. - * \returns true if the point is in the triangle, false if it is outside. + * \returns True if the point is in the triangle, false if it is outside. */ bool IsInTriangle(su2double val_CV1, su2double val_CV2, unsigned long val_id_triangle, unsigned long i_level = 0); /*! - * \brief compute the area of a triangle given the 3 points of the triangle. - * \param[in] x1 - the coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). - * \param[in] y1 - the coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). - * \param[in] x2 - the coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). - * \param[in] y2 - the coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). - * \param[in] x3 - the coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). - * \param[in] y3 - the coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). - * \returns the absolute value of the area of the triangle. + * \brief Compute the area of a triangle given the 3 points of the triangle. + * \param[in] x1 - The coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). + * \param[in] y1 - The coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). + * \param[in] x2 - The coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). + * \param[in] y2 - The coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). + * \param[in] x3 - The coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). + * \param[in] y3 - The coordinates of the points P1(x1,y1), P2(x2,y2) and P3(x3,y3). + * \returns The absolute value of the area of the triangle. */ inline su2double TriArea(su2double x1, su2double y1, su2double x2, su2double y2, su2double x3, su2double y3) { return abs((x1 * (y2 - y3) + x2 * (y3 - y1) + x3 * (y1 - y2)) * 0.5); @@ -261,10 +273,10 @@ class CLookUpTable { /*! * \brief Compute the values of the first and second controlling variable based on normalized query coordinates * \param[in] inclusion_levels - Pair containing lower(first) and upper(second) table inclusion level indices. - * \param[in] val_CV1 - first controlling variable value. - * \param[in] val_CV2 - second controlling variable value. - * \param[in] val_CV3 - third controlling variable value. - * \returns array with first and second controlling variable values for the upper and lower table inclusion levels. + * \param[in] val_CV1 - First controlling variable value. + * \param[in] val_CV2 - Second controlling variable value. + * \param[in] val_CV3 - Third controlling variable value. + * \returns Array with first and second controlling variable values for the upper and lower table inclusion levels. */ std::array, 2> ComputeNormalizedXY(std::pair& inclusion_levels, const su2double val_CV1, const su2double val_CV2, @@ -274,60 +286,63 @@ class CLookUpTable { CLookUpTable(const std::string& file_name_lut, std::string name_CV1_in, std::string name_CV2_in); /*! - * \brief print information to screen. + * \brief Print information to screen. */ void PrintTableInfo(); /*! - * \brief lookup 1 value of the single variable "val_name_var" using controlling variable values(val_CV1,val_CV2). - * \param[in] val_name_var - string name of the variable to look up. - * \param[out] val_var - the stored value of the variable to look up. - * \param[in] val_CV1 - value of controlling variable 1. - * \param[in] val_CV2 - value of controlling variable 2. + * \brief Lookup 1 value of the single variable "val_name_var" using controlling variable values(val_CV1,val_CV2). + * \param[in] val_name_var - String name of the variable to look up. + * \param[out] val_var - The stored value of the variable to look up. + * \param[in] val_CV1 - Value of controlling variable 1. + * \param[in] val_CV2 - Value of controlling variable 2. * \returns 1 if the lookup and subsequent interpolation was a success, 0 if not. */ unsigned long LookUp_XY(const std::string& val_name_var, su2double* val_var, su2double val_CV1, su2double val_CV2, unsigned long i_level = 0); /*! - * \brief lookup 1 value for each of the variables in "val_name_var" using controlling variable - * values(val_CV1,val_CV2). \param[in] val_names_var - vector of string names of the variables to look up. \param[out] - * val_vars - pointer to the vector of stored values of the variables to look up. \param[in] val_CV1 - value of - * controlling variable 1. \param[in] val_CV2 - value of controlling variable 2. \returns 1 if the lookup and - * subsequent interpolation was a success, 0 if not. + * \brief Lookup 1 value for each of the variables in "val_name_var" using controlling variable + * values(val_CV1,val_CV2). + * \param[in] val_names_var - vector of string names of the variables to look up. + * \param[out] val_vars - pointer to the vector of stored values of the variables to look up. + * \param[in] val_CV1 - value of controlling variable 1. + * \param[in] val_CV2 - value of controlling variable 2. + * \returns 1 if the lookup and subsequent interpolation was a success, 0 if not. */ unsigned long LookUp_XY(const std::vector& val_names_var, std::vector& val_vars, su2double val_CV1, su2double val_CV2, unsigned long i_level = 0); /*! - * \brief lookup the value of the variable "val_name_var" using controlling variable values(val_CV1,val_CV2). - * \param[in] val_name_var - string name of the variable to look up. - * \param[out] val_var - the stored value of the variable to look up. - * \param[in] val_CV1 - value of controlling variable 1. - * \param[in] val_CV2 - value of controlling variable 2. + * \brief Lookup the value of the variable "val_name_var" using controlling variable values(val_CV1,val_CV2). + * \param[in] val_name_var - String name of the variable to look up. + * \param[out] val_var - The stored value of the variable to look up. + * \param[in] val_CV1 - Value of controlling variable 1. + * \param[in] val_CV2 - Value of controlling variable 2. * \returns 1 if the lookup and subsequent interpolation was a success, 0 if not. */ unsigned long LookUp_XY(const std::vector& val_names_var, std::vector& val_vars, su2double val_CV1, su2double val_CV2, unsigned long i_level = 0); + /*! - * \brief lookup the value of the variable "val_name_var" using controlling variable values(val_CV1,val_CV2,val_z). - * \param[in] val_name_var - string name of the variable to look up. - * \param[out] val_var - the stored value of the variable to look up. - * \param[in] val_CV1 - value of controlling variable 1. - * \param[in] val_CV2 - value of controlling variable 2. - * \param[in] val_CV3 - value of controlling variable 3. + * \brief Lookup the value of the variable "val_name_var" using controlling variable values(val_CV1,val_CV2,val_z). + * \param[in] val_name_var - String name of the variable to look up. + * \param[out] val_var - The stored value of the variable to look up. + * \param[in] val_CV1 - Value of controlling variable 1. + * \param[in] val_CV2 - Value of controlling variable 2. + * \param[in] val_CV3 - Value of controlling variable 3. * \returns 1 if the lookup and subsequent interpolation was a success, 0 if not. */ unsigned long LookUp_XYZ(const std::string& val_name_var, su2double* val_var, su2double val_CV1, su2double val_CV2, su2double val_CV3); /*! - * \brief lookup the value of the variable "val_name_var" using controlling variable values(val_CV1,val_CV2,val_z). - * \param[in] val_name_var - string name of the variable to look up. - * \param[out] val_var - the stored value of the variable to look up. - * \param[in] val_CV1 - value of controlling variable 1. - * \param[in] val_CV2 - value of controlling variable 2. - * \param[in] val_CV3 - value of controlling variable 3. + * \brief Lookup the value of the variable "val_name_var" using controlling variable values(val_CV1,val_CV2,val_z). + * \param[in] val_name_var - String name of the variable to look up. + * \param[out] val_var - The stored value of the variable to look up. + * \param[in] val_CV1 - Value of controlling variable 1. + * \param[in] val_CV2 - Value of controlling variable 2. + * \param[in] val_CV3 - Value of controlling variable 3. * \returns 1 if the lookup and subsequent interpolation was a success, 0 if not. */ unsigned long LookUp_XYZ(const std::vector& val_names_var, std::vector& val_vars, @@ -335,23 +350,23 @@ class CLookUpTable { /*! * \brief Find the table levels with constant z-values directly above and below query val_z. - * \param[in] val_CV3 - value of controlling variable 3. + * \param[in] val_CV3 - Value of controlling variable 3. * \param[in] within_limits - Whether query point lies within table bounds. - * \returns pair of inclusion level indices (first = lower level index, second = upper level index) + * \returns Pair of inclusion level indices (first = lower level index, second = upper level index). */ std::pair FindInclusionLevels(const su2double val_CV3); /*! - * \brief determine the minimum and maximum value of the second controlling variable. - * \returns pair of minimum and maximum value of controlling variable 2. + * \brief Determine the minimum and maximum value of the second controlling variable. + * \returns Pair of minimum and maximum value of controlling variable 2. */ inline std::pair GetTableLimitsY(unsigned long i_level = 0) const { return limits_table_y[i_level]; } /*! - * \brief determine the minimum and maximum value of the first controlling variable. - * \returns pair of minimum and maximum value of controlling variable 1. + * \brief Determine the minimum and maximum value of the first controlling variable. + * \returns Pair of minimum and maximum value of controlling variable 1. */ inline std::pair GetTableLimitsX(unsigned long i_level = 0) const { return limits_table_x[i_level]; diff --git a/Common/include/option_structure.hpp b/Common/include/option_structure.hpp index d8a6e4dd69e..ea657b43afd 100644 --- a/Common/include/option_structure.hpp +++ b/Common/include/option_structure.hpp @@ -472,28 +472,24 @@ enum RUNTIME_TYPE { RUNTIME_ADJSPECIES_SYS = 26,/*!< \brief One-physics case, the code is solving the adjoint species model. */ }; -const int FLOW_SOL = 0; /*!< \brief Position of the mean flow solution in the solver container array. */ -const int ADJFLOW_SOL = 1; /*!< \brief Position of the continuous adjoint flow solution in the solver container array. */ - -const int TURB_SOL = 2; /*!< \brief Position of the turbulence model solution in the solver container array. */ -const int ADJTURB_SOL = 3; /*!< \brief Position of the continuous adjoint turbulence solution in the solver container array. */ - -const int TRANS_SOL = 4; /*!< \brief Position of the transition model solution in the solver container array. */ -const int HEAT_SOL = 5; /*!< \brief Position of the heat equation in the solution solver array. */ -const int ADJHEAT_SOL = 6; /*!< \brief Position of the adjoint heat equation in the solution solver array. */ -const int RAD_SOL = 7; /*!< \brief Position of the radiation equation in the solution solver array. */ -const int ADJRAD_SOL = 8; /*!< \brief Position of the continuous adjoint turbulence solution in the solver container array. */ - -const int MESH_SOL = 9; /*!< \brief Position of the mesh solver. */ -const int ADJMESH_SOL = 10; /*!< \brief Position of the adjoint of the mesh solver. */ - -const int SPECIES_SOL = 11; /*!< \brief Position of the species solver. */ -const int ADJSPECIES_SOL = 12; /*!< \brief Position of the adjoint of the species solver. */ - -const int FEA_SOL = 0; /*!< \brief Position of the FEA equation in the solution solver array. */ -const int ADJFEA_SOL = 1; /*!< \brief Position of the FEA adjoint equation in the solution solver array. */ - -const int TEMPLATE_SOL = 0; /*!< \brief Position of the template solution. */ + enum SOLVER_TYPE : const int { + FLOW_SOL=0, /*!< \brief Position of the mean flow solution in the solver container array. */ + ADJFLOW_SOL=1, /*!< \brief Position of the continuous adjoint flow solution in the solver container array. */ + TURB_SOL=2, /*!< \brief Position of the turbulence model solution in the solver container array. */ + ADJTURB_SOL=3, /*!< \brief Position of the continuous adjoint turbulence solution in the solver container array. */ + TRANS_SOL=4, /*!< \brief Position of the transition model solution in the solver container array. */ + HEAT_SOL=5, /*!< \brief Position of the heat equation in the solution solver array. */ + ADJHEAT_SOL=6, /*!< \brief Position of the adjoint heat equation in the solution solver array. */ + RAD_SOL=7, /*!< \brief Position of the radiation equation in the solution solver array. */ + ADJRAD_SOL=8, /*!< \brief Position of the continuous adjoint turbulence solution in the solver container array. */ + MESH_SOL=9, /*!< \brief Position of the mesh solver. */ + ADJMESH_SOL=10, /*!< \brief Position of the adjoint of the mesh solver. */ + SPECIES_SOL=11, /*!< \brief Position of the species solver. */ + ADJSPECIES_SOL=12, /*!< \brief Position of the adjoint of the species solver. */ + FEA_SOL=0, /*!< \brief Position of the Finite Element flow solution in the solver container array. */ + ADJFEA_SOL=1, /*!< \brief Position of the continuous adjoint Finite Element flow solution in the solver container array. */ + TEMPLATE_SOL=0, /*!< \brief Position of the template solution. */ + }; const int CONV_TERM = 0; /*!< \brief Position of the convective terms in the numerics container array. */ const int VISC_TERM = 1; /*!< \brief Position of the viscous terms in the numerics container array. */ @@ -551,6 +547,7 @@ enum ENUM_FLUIDMODEL { SU2_NONEQ = 8, /*!< \brief User defined gas model for nonequilibrium flow. */ FLUID_MIXTURE = 9, /*!< \brief Species mixture model. */ COOLPROP = 10, /*!< \brief Thermodynamics library. */ + FLUID_FLAMELET = 11, /*!< \brief lookup table (LUT) method for premixed flamelets. */ }; static const MapType FluidModel_Map = { MakePair("STANDARD_AIR", STANDARD_AIR) @@ -564,6 +561,7 @@ static const MapType FluidModel_Map = { MakePair("SU2_NONEQ", SU2_NONEQ) MakePair("FLUID_MIXTURE", FLUID_MIXTURE) MakePair("COOLPROP", COOLPROP) + MakePair("FLUID_FLAMELET", FLUID_FLAMELET) }; /*! @@ -653,12 +651,14 @@ enum class VISCOSITYMODEL { CONSTANT, /*!< \brief Constant viscosity. */ SUTHERLAND, /*!< \brief Sutherlands Law viscosity. */ POLYNOMIAL, /*!< \brief Polynomial viscosity. */ + FLAMELET, /*!< \brief LUT method for flamelets */ COOLPROP, /*!< \brief CoolProp viscosity. */ }; static const MapType ViscosityModel_Map = { MakePair("CONSTANT_VISCOSITY", VISCOSITYMODEL::CONSTANT) MakePair("SUTHERLAND", VISCOSITYMODEL::SUTHERLAND) MakePair("POLYNOMIAL_VISCOSITY", VISCOSITYMODEL::POLYNOMIAL) + MakePair("FLAMELET", VISCOSITYMODEL::FLAMELET) MakePair("COOLPROP", VISCOSITYMODEL::COOLPROP) }; @@ -681,12 +681,14 @@ enum class CONDUCTIVITYMODEL { CONSTANT, /*!< \brief Constant thermal conductivity. */ CONSTANT_PRANDTL, /*!< \brief Constant Prandtl number. */ POLYNOMIAL, /*!< \brief Polynomial thermal conductivity. */ + FLAMELET, /*!< \brief LUT method for flamelets */ COOLPROP, /*!< \brief COOLPROP thermal conductivity. */ }; static const MapType ConductivityModel_Map = { MakePair("CONSTANT_CONDUCTIVITY", CONDUCTIVITYMODEL::CONSTANT) MakePair("CONSTANT_PRANDTL", CONDUCTIVITYMODEL::CONSTANT_PRANDTL) MakePair("POLYNOMIAL_CONDUCTIVITY", CONDUCTIVITYMODEL::POLYNOMIAL) + MakePair("FLAMELET", CONDUCTIVITYMODEL::FLAMELET) MakePair("COOLPROP", CONDUCTIVITYMODEL::COOLPROP) }; @@ -710,6 +712,7 @@ enum class DIFFUSIVITYMODEL { CONSTANT_SCHMIDT, /*!< \brief Constant Schmidt number for mass diffusion in scalar transport. */ UNITY_LEWIS, /*!< \brief Unity Lewis model for mass diffusion in scalar transport. */ CONSTANT_LEWIS, /*!< \brief Different Lewis number model for mass diffusion in scalar transport. */ + FLAMELET, /*!< \brief flamelet model for tabulated chemistry, diffusivity from lookup table */ }; static const MapType Diffusivity_Model_Map = { @@ -717,6 +720,7 @@ static const MapType Diffusivity_Model_Map = { MakePair("CONSTANT_SCHMIDT", DIFFUSIVITYMODEL::CONSTANT_SCHMIDT) MakePair("UNITY_LEWIS", DIFFUSIVITYMODEL::UNITY_LEWIS) MakePair("CONSTANT_LEWIS", DIFFUSIVITYMODEL::CONSTANT_LEWIS) + MakePair("FLAMELET", DIFFUSIVITYMODEL::FLAMELET) }; /*! @@ -1287,11 +1291,29 @@ inline LM_ParsedOptions ParseLMOptions(const LM_OPTIONS *LM_Options, unsigned sh */ enum class SPECIES_MODEL { NONE, /*!< \brief No scalar transport model. */ - SPECIES_TRANSPORT, /*!< \brief Passive scalar transport model. */ + SPECIES_TRANSPORT, /*!< \brief species transport model. */ + FLAMELET, /*!< \brief flamelet model. */ }; static const MapType Species_Model_Map = { MakePair("NONE", SPECIES_MODEL::NONE) MakePair("SPECIES_TRANSPORT", SPECIES_MODEL::SPECIES_TRANSPORT) + MakePair("FLAMELET", SPECIES_MODEL::FLAMELET) +}; + +/*! + * \brief Progress variable and enthalpy are the first and second entries in the lookup table. + * \note The order matters. + */ +enum FLAMELET_SCALAR_VARIABLES { + I_PROGVAR, + I_ENTH, +}; + +/*! + * \brief the source terms for the flamelet method. At the moment only progress variable + */ +enum FLAMELET_SCALAR_SOURCES { + I_SRC_TOT_PROGVAR }; /*! diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index e705b7f38d4..ace96e28cbb 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -837,6 +837,7 @@ void CConfig::SetPointersNull() { Marker_Designing = nullptr; Marker_GeoEval = nullptr; Marker_Plotting = nullptr; Marker_Analyze = nullptr; Marker_PyCustom = nullptr; Marker_WallFunctions = nullptr; Marker_CfgFile_KindBC = nullptr; Marker_All_KindBC = nullptr; Marker_SobolevBC = nullptr; + Marker_StrongBC = nullptr; Kind_WallFunctions = nullptr; IntInfo_WallFunctions = nullptr; @@ -1351,8 +1352,15 @@ void CConfig::SetConfig_Options() { addDoubleListOption("SPECIES_CLIPPING_MAX", nSpecies_Clipping_Max, Species_Clipping_Max); /*!\brief SPECIES_CLIPPING_MIN \n DESCRIPTION: Minimum values for scalar clipping \ingroup Config*/ addDoubleListOption("SPECIES_CLIPPING_MIN", nSpecies_Clipping_Min, Species_Clipping_Min); - /*!\brief SPECIES_CLIPPING \n DESCRIPTION: Use strong inlet and outlet BC in the species solver \n DEFAULT: false \ingroup Config*/ - addBoolOption("SPECIES_USE_STRONG_BC", Species_StrongBC, false); + + /*!\brief FLAME_INIT \n DESCRIPTION: flame initialization using the flamelet model \ingroup Config*/ + /*--- flame offset (x,y,z) ---*/ + flame_init[0] = 0.0; flame_init[1] = 0.0; flame_init[2] = 0.0; + /*--- flame normal (nx, ny, nz) ---*/ + flame_init[3] = 1.0; flame_init[4] = 0.0; flame_init[5] = 0.0; + /*--- flame thickness (x) and flame burnt thickness (after this thickness, we have unburnt conditions again) ---*/ + flame_init[6] = 0.5e-3; flame_init[7] = 1.0; + addDoubleArrayOption("FLAME_INIT", 8,flame_init); /*--- Options related to mass diffusivity and thereby the species solver. ---*/ @@ -1496,6 +1504,9 @@ void CConfig::SetConfig_Options() { addWallFunctionOption("MARKER_WALL_FUNCTIONS", nMarker_WallFunctions, Marker_WallFunctions, Kind_WallFunctions, IntInfo_WallFunctions, DoubleInfo_WallFunctions); + /*!\brief MARKER_STRONG_BC\n DESCRIPTION: Markers where a strong BC must be applied.*/ + addStringListOption("MARKER_SPECIES_STRONG_BC", nMarker_StrongBC, Marker_StrongBC); + /*!\brief ACTDISK_TYPE \n DESCRIPTION: Actuator Disk boundary type \n OPTIONS: see \link ActDisk_Map \endlink \n Default: VARIABLES_JUMP \ingroup Config*/ addEnumOption("ACTDISK_TYPE", Kind_ActDisk, ActDisk_Map, VARIABLES_JUMP); @@ -2082,6 +2093,18 @@ void CConfig::SetConfig_Options() { /* DESCRIPTION: Determine if we need to allocate memory to store the multizone residual. \n DEFAULT: true (temporarily) */ addBoolOption("MULTIZONE_RESIDUAL", Multizone_Residual, false); + /*!\brief File name of the flamelet look up table.*/ + addStringOption("FILENAME_LUT", file_name_lut, string("LUT")); + + /* DESCRIPTION: Names of the passive lookup variables for flamelet LUT */ + addStringListOption("LOOKUP_NAMES", n_lookups, table_lookup_names); + + /* DESCRIPTION: Names of the user transport equations solved in the flamelet problem. */ + addStringListOption("USER_SCALAR_NAMES", n_user_scalars, user_scalar_names); + + /* DESCRIPTION: Names of the user scalar source terms. */ + addStringListOption("USER_SOURCE_NAMES", n_user_sources, user_source_names); + /*!\brief CONV_FILENAME \n DESCRIPTION: Output file convergence history (w/o extension) \n DEFAULT: history \ingroup Config*/ addStringOption("CONV_FILENAME", Conv_FileName, string("history")); /*!\brief BREAKDOWN_FILENAME \n DESCRIPTION: Output file forces breakdown \ingroup Config*/ @@ -2987,6 +3010,8 @@ void CConfig::SetConfig_Parsing(istream& config_buffer){ newString.append("DYN_TIME is deprecated. Use MAX_TIME instead.\n\n"); else if (!option_name.compare("DYNAMIC_ANALYSIS")) newString.append("DYNAMIC_ANALYSIS is deprecated. Use TIME_DOMAIN instead.\n\n"); + else if (!option_name.compare("SPECIES_USE_STRONG_BC")) + newString.append("SPECIES_USE_STRONG_BC is deprecated. Use MARKER_SPECIES_STRONG_BC= (marker1, ...) instead.\n\n"); else { /*--- Find the most likely candidate for the unrecognized option, based on the length of start and end character sequences shared by candidates and the option. ---*/ @@ -3297,6 +3322,7 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i (Kind_FluidModel == IDEAL_GAS) || (Kind_FluidModel == INC_IDEAL_GAS) || (Kind_FluidModel == FLUID_MIXTURE) || + (Kind_FluidModel == FLUID_FLAMELET) || (Kind_FluidModel == INC_IDEAL_GAS_POLY) || (Kind_FluidModel == CONSTANT_DENSITY)); bool noneq_gas = ((Kind_FluidModel == MUTATIONPP) || @@ -3644,9 +3670,11 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i } } - if (Kind_ObjFunc[0] == CUSTOM_OBJFUNC && CustomObjFunc.empty() && !Multizone_Problem) { - SU2_MPI::Error("The expression for the custom objective function was not set.\n" - "For example, CUSTOM_OBJFUNC= LIFT/DRAG", CURRENT_FUNCTION); + if (nObj > 0){ + if (Kind_ObjFunc[0] == CUSTOM_OBJFUNC && CustomObjFunc.empty() && !Multizone_Problem) { + SU2_MPI::Error("The expression for the custom objective function was not set.\n" + "For example, CUSTOM_OBJFUNC= LIFT/DRAG", CURRENT_FUNCTION); + } } /*--- Check for unsteady problem ---*/ @@ -3845,7 +3873,7 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i } break; default: - if (nSpecies_Init + 1 != 1) SU2_MPI::Error("Viscosity model not available.", CURRENT_FUNCTION); + if (nSpecies_Init + 1 != 1) SU2_MPI::Error("Fluid mixture: viscosity model not available.", CURRENT_FUNCTION); break; } @@ -3892,6 +3920,35 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i } } + + if (Kind_Species_Model == SPECIES_MODEL::FLAMELET) { + + if (Kind_FluidModel != FLUID_FLAMELET) { + SU2_MPI::Error("The use of SCALAR_MODEL= FLAMELET requires the FLUID_MODEL option to be FLUID_FLAMELET", + CURRENT_FUNCTION); + } + + if (Kind_DensityModel != INC_DENSITYMODEL::VARIABLE) { + SU2_MPI::Error("The use of FLUID_FLAMELET requires the INC_DENSITY_MODEL option to be VARIABLE", + CURRENT_FUNCTION); + } + + if (Kind_ConductivityModel != CONDUCTIVITYMODEL::FLAMELET) { + SU2_MPI::Error("The use of FLUID_FLAMELET requires the CONDUCTIVITY_MODEL option to be FLAMELET", + CURRENT_FUNCTION); + } + + if (Kind_Diffusivity_Model != DIFFUSIVITYMODEL::FLAMELET) { + SU2_MPI::Error("The use of FLUID_FLAMELET requires the DIFFUSIVITY_MODEL option to be FLAMELET", + CURRENT_FUNCTION); + } + + if (Kind_ViscosityModel != VISCOSITYMODEL::FLAMELET) { + SU2_MPI::Error("The use of FLUID_FLAMELET requires the VISCOSITY_MODEL option to be FLAMELET", + CURRENT_FUNCTION); + } + } + /*--- Check for Measurement System ---*/ if (SystemMeasurements == US && !standard_air) { @@ -4859,9 +4916,9 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i Kind_FluidModel = CONSTANT_DENSITY; } - /*--- Energy equation must be active for any fluid models other than constant density. ---*/ - if ((Kind_DensityModel != INC_DENSITYMODEL::CONSTANT) && (Kind_Species_Model==SPECIES_MODEL::NONE)) Energy_Equation = true; + /*--- For the flamelet combustion model, energy equation is a passive field, we lookup T and write it to the field ---*/ + if (Kind_Species_Model == SPECIES_MODEL::FLAMELET ) Energy_Equation = false; if (Kind_DensityModel == INC_DENSITYMODEL::BOUSSINESQ) { Energy_Equation = true; @@ -4871,7 +4928,7 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i } if (Kind_DensityModel == INC_DENSITYMODEL::VARIABLE) { - if (Kind_FluidModel != INC_IDEAL_GAS && Kind_FluidModel != INC_IDEAL_GAS_POLY && Kind_FluidModel != FLUID_MIXTURE) { + if (Kind_FluidModel != INC_IDEAL_GAS && Kind_FluidModel != INC_IDEAL_GAS_POLY && Kind_FluidModel != FLUID_MIXTURE && Kind_FluidModel != FLUID_FLAMELET) { SU2_MPI::Error("Variable density incompressible solver limited to ideal gases.\n Check the fluid model options (use INC_IDEAL_GAS, INC_IDEAL_GAS_POLY).", CURRENT_FUNCTION); } } @@ -5313,7 +5370,7 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i CURRENT_FUNCTION); /*--- Checks for additional species transport. ---*/ - if (Kind_Species_Model == SPECIES_MODEL::SPECIES_TRANSPORT) { + if ((Kind_Species_Model == SPECIES_MODEL::SPECIES_TRANSPORT) || (Kind_Species_Model == SPECIES_MODEL::FLAMELET)) { if (Kind_Solver != MAIN_SOLVER::INC_NAVIER_STOKES && Kind_Solver != MAIN_SOLVER::INC_RANS && Kind_Solver != MAIN_SOLVER::DISC_ADJ_INC_NAVIER_STOKES && @@ -5383,24 +5440,35 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i nSpecies = nSpecies_Init; /*--- Check whether some variables (or their sums) are in physical bounds. [0,1] for species related quantities. ---*/ - su2double Species_Init_Sum = 0.0; - for (unsigned short iSpecies = 0; iSpecies < nSpecies; iSpecies++) { - checkScalarBounds(Species_Init[iSpecies], "SPECIES_INIT individual", 0.0, 1.0); - Species_Init_Sum += Species_Init[iSpecies]; - } - checkScalarBounds(Species_Init_Sum, "SPECIES_INIT sum", 0.0, 1.0); - - for (unsigned short iMarker = 0; iMarker < nMarker_Inlet_Species; iMarker++) { - su2double Inlet_SpeciesVal_Sum = 0.0; + /*--- 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++) { - checkScalarBounds(Inlet_SpeciesVal[iMarker][iSpecies], "MARKER_INLET_SPECIES individual", 0.0, 1.0); - Inlet_SpeciesVal_Sum += Inlet_SpeciesVal[iMarker][iSpecies]; + checkScalarBounds(Species_Init[iSpecies], "SPECIES_INIT individual", 0.0, 1.0); + Species_Init_Sum += Species_Init[iSpecies]; + } + checkScalarBounds(Species_Init_Sum, "SPECIES_INIT sum", 0.0, 1.0); + + for (iMarker = 0; iMarker < nMarker_Inlet_Species; iMarker++) { + su2double Inlet_SpeciesVal_Sum = 0.0; + for (unsigned short iSpecies = 0; iSpecies < nSpecies; iSpecies++) { + checkScalarBounds(Inlet_SpeciesVal[iMarker][iSpecies], "MARKER_INLET_SPECIES individual", 0.0, 1.0); + Inlet_SpeciesVal_Sum += Inlet_SpeciesVal[iMarker][iSpecies]; + } + checkScalarBounds(Inlet_SpeciesVal_Sum, "MARKER_INLET_SPECIES sum", 0.0, 1.0); } - checkScalarBounds(Inlet_SpeciesVal_Sum, "MARKER_INLET_SPECIES sum", 0.0, 1.0); } } // species transport checks + /*--- Define some variables for flamelet model. ---*/ + if (Kind_Species_Model == SPECIES_MODEL::FLAMELET) { + /*--- The controlling variables are progress variable and total enthalpy ---*/ + n_control_vars = 2; + /*--- We can have additional user defined transported scalars ---*/ + n_scalars = n_control_vars + n_user_scalars; + } + if (Kind_Regime == ENUM_REGIME::COMPRESSIBLE && GetBounded_Scalar()) { SU2_MPI::Error("BOUNDED_SCALAR discretization can only be used for incompressible problems.", CURRENT_FUNCTION); } @@ -5946,7 +6014,7 @@ void CConfig::SetOutput(SU2_COMPONENT val_software, unsigned short val_izone) { iMarker_Designing, iMarker_GeoEval, iMarker_Plotting, iMarker_Analyze, iMarker_DV, iDV_Value, iMarker_ZoneInterface, iMarker_PyCustom, iMarker_Load_Dir, iMarker_Disp_Dir, iMarker_Load_Sine, iMarker_Clamped, iMarker_Moving, iMarker_Supersonic_Inlet, iMarker_Supersonic_Outlet, iMarker_ActDiskInlet, - iMarker_Emissivity, + iMarker_Emissivity, iMarker_StrongBC, iMarker_ActDiskOutlet, iMarker_MixingPlaneInterface, iMarker_SobolevBC; @@ -7458,6 +7526,15 @@ void CConfig::SetOutput(SU2_COMPONENT val_software, unsigned short val_izone) { BoundaryTable.PrintFooter(); } + if (nMarker_StrongBC != 0) { + BoundaryTable << "Strong boundary"; + for (iMarker_StrongBC = 0; iMarker_StrongBC < nMarker_StrongBC; iMarker_StrongBC++) { + BoundaryTable << Marker_StrongBC[iMarker_StrongBC]; + if (iMarker_StrongBC < nMarker_StrongBC-1) BoundaryTable << " "; + } + BoundaryTable.PrintFooter(); + } + if (nMarker_Custom != 0) { BoundaryTable << "Custom boundary"; for (iMarker_Custom = 0; iMarker_Custom < nMarker_Custom; iMarker_Custom++) { @@ -9493,6 +9570,15 @@ su2double CConfig::GetWall_Emissivity(const string& val_marker) const { return Wall_Emissivity[iMarker_Emissivity]; } +bool CConfig::GetMarker_StrongBC(const string& val_marker) const { + + unsigned short iMarker_StrongBC = 0; + + for (iMarker_StrongBC = 0; iMarker_StrongBC < nMarker_StrongBC; iMarker_StrongBC++) + if (Marker_StrongBC[iMarker_StrongBC] == val_marker) return true; + + return false; +} su2double CConfig::GetFlowLoad_Value(const string& val_marker) const { unsigned short iMarker_FlowLoad; for (iMarker_FlowLoad = 0; iMarker_FlowLoad < nMarker_FlowLoad; iMarker_FlowLoad++) diff --git a/Common/src/containers/CFileReaderLUT.cpp b/Common/src/containers/CFileReaderLUT.cpp index d71f4e9ee3a..31bb8de3338 100644 --- a/Common/src/containers/CFileReaderLUT.cpp +++ b/Common/src/containers/CFileReaderLUT.cpp @@ -1,5 +1,5 @@ /*! - * \file CFileReaderLUT.hpp + * \file CFileReaderLUT.cpp * \brief reading lookup table for tabulated fluid properties * \author D. Mayer, T. Economon * \version 7.5.1 "Blackbird" diff --git a/Common/src/containers/CLookUpTable.cpp b/Common/src/containers/CLookUpTable.cpp index 130f53bfc18..d13feb6b82c 100644 --- a/Common/src/containers/CLookUpTable.cpp +++ b/Common/src/containers/CLookUpTable.cpp @@ -92,15 +92,14 @@ CLookUpTable::CLookUpTable(const string& var_file_name_lut, string name_CV1_in, if (rank == MASTER_NODE) { switch (table_dim) { case 2: - cout << "Construction of trapezoidal map took " << stopTime - startTime << " seconds\n" << endl; + cout << "\nConstruction of trapezoidal map took " << stopTime - startTime << " seconds\n" << endl; break; case 3: - cout << "Construction of trapezoidal map stack took " << stopTime - startTime << " seconds\n" << endl; + cout << "\nConstruction of trapezoidal map stack took " << stopTime - startTime << " seconds\n" << endl; break; default: break; } - cout << "Precomputing interpolation coefficients..." << endl; } ComputeInterpCoeffs(); @@ -576,7 +575,7 @@ unsigned long CLookUpTable::LookUp_XY(const string& val_name_var, su2double* val su2double val_CV2, unsigned long i_level) { unsigned long exit_code = 1; - if (val_name_var.compare("NULL") == 0) { + if (noSource(val_name_var)) { *val_var = 0.0; exit_code = 0; return exit_code; @@ -645,18 +644,19 @@ unsigned long CLookUpTable::LookUp_XY(const vector& val_names_var, vecto } /* loop over variable names and interpolate / get values */ - if (exit_code == 0) { - for (long unsigned int i_var = 0; i_var < val_names_var.size(); ++i_var) { - if (val_names_var[i_var].compare("NULL") == 0) { - *val_vars[i_var] = 0.0; - } else { + for (long unsigned int i_var = 0; i_var < val_names_var.size(); ++i_var) { + if (noSource(val_names_var[i_var])) { + *val_vars[i_var] = 0.0; + exit_code = 0; + } else { + if (exit_code == 0) { /* first, copy the single triangle from the large triangle list*/ for (int p = 0; p < 3; p++) triangle[p] = triangles[i_level][id_triangle][p]; *val_vars[i_var] = Interpolate(GetDataP(val_names_var[i_var], i_level), triangle, interp_coeffs); + } else { + InterpolateToNearestNeighbors(val_CV1, val_CV2, val_names_var[i_var], val_vars[i_var], i_level); } } - } else { - InterpolateToNearestNeighbors(val_CV1, val_CV2, val_names_var, val_vars, i_level); } return exit_code; diff --git a/Common/src/containers/CTrapezoidalMap.cpp b/Common/src/containers/CTrapezoidalMap.cpp index 22f45868a8f..c0f71e98eb0 100644 --- a/Common/src/containers/CTrapezoidalMap.cpp +++ b/Common/src/containers/CTrapezoidalMap.cpp @@ -26,6 +26,7 @@ */ #include +#include #include "../../Common/include/option_structure.hpp" #include "../../Common/include/containers/CTrapezoidalMap.hpp" @@ -40,6 +41,9 @@ using namespace std; CTrapezoidalMap::CTrapezoidalMap(const su2double* samples_x, const su2double* samples_y, const unsigned long size, vector > const& edges, su2vector > const& val_edge_to_triangle) { + int rank = SU2_MPI::GetRank(); + su2double startTime = SU2_MPI::Wtime(); + edge_to_triangle = su2vector >(val_edge_to_triangle); unique_bands_x.assign(samples_x, samples_x + size); @@ -71,6 +75,7 @@ CTrapezoidalMap::CTrapezoidalMap(const su2double* samples_x, const su2double* sa unsigned long n_edges = edges.size(); /* edge index */ unsigned long i_edge = 0; + unsigned long j_edge = 0; /* counter for edges intersects */ unsigned long n_intersects = 0; /* lower and upper x value of each band */ @@ -125,6 +130,52 @@ CTrapezoidalMap::CTrapezoidalMap(const su2double* samples_x, const su2double* sa i_band++; } + + su2double stopTime = SU2_MPI::Wtime(); + + /* calculate size of trapezoidal map components */ + double size_unique_bands = sizeof(su2double) * unique_bands_x.size() / 1e6; + double size_edge_limits_x = sizeof(su2double) * edge_limits_x.size() * 2 / 1e6; + double size_edge_limits_y = sizeof(su2double) * edge_limits_y.size() * 2 / 1e6; + + double size_edge_to_triangle = 0; + for (i_edge = 0; i_edge < edge_to_triangle.size(); i_edge++) + for (j_edge = 0; j_edge < edge_to_triangle[i_edge].size(); j_edge++) + size_edge_to_triangle += sizeof(unsigned long) / 1e6; + + double size_y_edge_at_band_mid = 0; + for (unsigned long i_y = 0; i_y < y_edge_at_band_mid.size(); i_y++) + for (unsigned long j_y = 0; j_y < y_edge_at_band_mid[i_y].size(); j_y++) + size_y_edge_at_band_mid += sizeof(su2double) / 1e6 + sizeof(unsigned long) / 1e6; + + double size_total = + size_unique_bands + size_edge_limits_x + size_edge_limits_y + size_edge_to_triangle + size_y_edge_at_band_mid; + + /* print size of trapezoidal map components to screen */ + if (rank == MASTER_NODE) { + cout << setfill(' '); + cout << "\n" << endl; + cout << "+------------------------------------------------------------------+\n"; + cout << "| Trapezoidal map info |\n"; + cout << "+------------------------------------------------------------------+" << endl; + + cout << "| Time to construct trapezoidal map: " << setw(22) << right << stopTime - startTime << " sec" + << " |" << endl; + cout << "| Size of unique_bands in memory: " << setw(22) << size_unique_bands << " MB " + << " |" << endl; + cout << "| Size of edge_limits_x in memory: " << setw(22) << size_edge_limits_x << " MB " + << " |" << endl; + cout << "| Size of edge_limits_y in memory: " << setw(22) << size_edge_limits_y << " MB " + << " |" << endl; + cout << "| Size of edge_to_triangle in memory: " << setw(22) << size_edge_to_triangle << " MB " + << " |" << endl; + cout << "| Size of y_edge_at_band_mid in memory: " << setw(22) << size_y_edge_at_band_mid << " MB " + << " |" << endl; + cout << "| Total: " << setw(22) << size_total << " MB " + << " |" << endl; + cout << "+------------------------------------------------------------------+" << endl; + cout << "\n" << endl; + } } unsigned long CTrapezoidalMap::GetTriangle(su2double val_x, su2double val_y) { @@ -134,19 +185,19 @@ unsigned long CTrapezoidalMap::GetTriangle(su2double val_x, su2double val_y) { /* within that band, find edges which enclose the (val_x, val_y) point */ pair edges = GetEdges(band, val_x, val_y); - /* identify the triangle using the two edges */ + /* identify the adjacent triangles using the two edges */ std::array triangles_edge_low; - - for (int i = 0; i < 2; i++) triangles_edge_low[i] = edge_to_triangle[edges.first][i]; + for (unsigned long i = 0; i < edge_to_triangle[edges.first].size(); i++) + triangles_edge_low[i] = edge_to_triangle[edges.first][i]; std::array triangles_edge_up; - for (int i = 0; i < 2; i++) triangles_edge_up[i] = edge_to_triangle[edges.second][i]; + for (unsigned long i = 0; i < edge_to_triangle[edges.second].size(); i++) + triangles_edge_up[i] = edge_to_triangle[edges.second][i]; sort(triangles_edge_low.begin(), triangles_edge_low.end()); sort(triangles_edge_up.begin(), triangles_edge_up.end()); - // The intersection of the faces to which upper or lower belongs is - // the face that both belong to. + /* The intersection of the faces to which upper or lower belongs is the face that both belong to. */ vector triangle(1); set_intersection(triangles_edge_up.begin(), triangles_edge_up.end(), triangles_edge_low.begin(), triangles_edge_low.end(), triangle.begin()); diff --git a/Common/src/grid_movement/CSurfaceMovement.cpp b/Common/src/grid_movement/CSurfaceMovement.cpp index 2147af42f9f..1905b57deaa 100644 --- a/Common/src/grid_movement/CSurfaceMovement.cpp +++ b/Common/src/grid_movement/CSurfaceMovement.cpp @@ -893,7 +893,7 @@ void CSurfaceMovement::SetParametricCoord(CGeometry* geometry, CConfig* config, cout << "Computed parametric coord for FFD box '" << FFDBox->GetTag() << "'\n"; cout << " Number of vertices (Total, Inside FFD, Mapped to FFD): " << GlobalVertex; cout << ", " << GlobalVisited << ", " << GlobalMapped << "\n"; - cout << " Max coord differece: " << MaxDiff << "\n"; + cout << " Max coord difference: " << MaxDiff << "\n"; } /*--- After the point inversion, copy the original information back. ---*/ diff --git a/SU2_CFD/include/fluid/CFluidFlamelet.hpp b/SU2_CFD/include/fluid/CFluidFlamelet.hpp new file mode 100644 index 00000000000..1d0bf91986d --- /dev/null +++ b/SU2_CFD/include/fluid/CFluidFlamelet.hpp @@ -0,0 +1,114 @@ +/*! + * \file CFluidFlamelet.hpp + * \brief Defines the flamelet fluid model + * \author D. Mayer, T. Economon, N. Beishuizen + * \version 7.5.1 "Blackbird" + * + * SU2 Project Website: https://su2code.github.io + * + * The SU2 Project is maintained by the SU2 Foundation + * (http://su2foundation.org) + * + * Copyright 2012-2023, SU2 Contributors (cf. AUTHORS.md) + * + * SU2 is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public + * License as published by the Free Software Foundation; either + * version 2.1 of the License, or (at your option) any later version. + * + * SU2 is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with SU2. If not, see . + */ + +#pragma once + +#include "../../Common/include/containers/CLookUpTable.hpp" +#include "CFluidModel.hpp" + +class CFluidFlamelet final : public CFluidModel { + private: + enum LOOKUP_TD { + TEMPERATURE, + DENSITY, + HEATCAPACITY, + VISCOSITY, + CONDUCTIVITY, + DIFFUSIONCOEFFICIENT, + MOLARWEIGHT, + SIZE + }; + + int rank; + + unsigned short n_scalars; + unsigned short n_lookups; + unsigned short n_user_scalars; /*!< \brief number of passive reactant species. */ + unsigned short n_control_vars; /*!< \brief number of controlling variables. */ + + vector table_scalar_names; /*!< \brief vector to store names of scalar variables. */ + vector table_lookup_names; /*!< \brief vector to store names of look up variables. */ + + su2double mass_diffusivity; /*!< \brief local mass diffusivity of the mixture */ + su2double molar_weight; /*!< \brief local molar weight of the mixture */ + + CLookUpTable* look_up_table; + + vector varnames_TD; /*!< \brief Lookup names for thermodynamic state variables. */ + vector val_vars_TD; /*!< \brief References to thermodynamic state variables. */ + + public: + CFluidFlamelet(CConfig* config, su2double value_pressure_operating); + + ~CFluidFlamelet(); + + /*! + * \brief Set the thermodynamic state. + * \param[in] val_temperature - temperature + * \param[in] val_scalars - pointer to species mass fractions + */ + void SetTDState_T(su2double val_temperature, const su2double* val_scalars = nullptr) override; + + /*! + * \brief Get the total enthalpy from the tabulated temperature and species (inverse lookup). + * \param[in/out] enthalpy - total enthalpy + * \param[in] val_prog - progress variable + * \param[in] val_temp - temperature + * \param[in] initial_value - initial value for the iterative lookup procedure. + * \param[out] exit_code = error code + */ + unsigned long GetEnthFromTemp(su2double& enthalpy, const su2double val_prog, const su2double val_temp, + su2double initial_value = 0) override; + + /*! + * \brief Return a pointer to the lookup table. + * \param[out] look_up_table - pointer to lookup table + */ + inline CLookUpTable* GetLookUpTable() override { return look_up_table; } + + /*! + * \brief Get the mass diffusivity of the species. + * \param[in] iVar - index to the species + * \param[out] mass_diffusivity - value of the mass diffusivity + */ + inline su2double GetMassDiffusivity(int iVar) override { return mass_diffusivity; } + + /*! + * \brief Get the thermal conductivity of the species. + * \param[in] iVar - index to the species + * \param[out] Kt - value of the thermal conductivity + */ + inline su2double GetThermalConductivity() override { return Kt; } + + /*! + * \brief Get the laminar viscosity of the species. + * \param[in] iVar - index to the species + * \param[out] Mu - value of the laminar viscosity + */ + inline su2double GetLaminarViscosity() override { return Mu; } + +}; diff --git a/SU2_CFD/include/fluid/CFluidModel.hpp b/SU2_CFD/include/fluid/CFluidModel.hpp index bd8b8b4952e..5e3785bcda1 100644 --- a/SU2_CFD/include/fluid/CFluidModel.hpp +++ b/SU2_CFD/include/fluid/CFluidModel.hpp @@ -38,6 +38,7 @@ using namespace std; +class CLookUpTable; /*! * \class CFluidModel * \brief Main class for defining the Thermo-Physical Model @@ -137,6 +138,24 @@ class CFluidModel { */ su2double GetCv() const { return Cv; } + /*! + * \brief Flamelet LUT - Get the number of transported scalars. + */ + virtual inline unsigned short GetNScalars() const { return 0; } + + /*! + * \brief Flamelet LUT - Get the lookup table. + */ + virtual CLookUpTable* GetLookUpTable() { return nullptr; } + + /*! + * \brief Flamelet LUT - Get the total enthalpy from the temperature (reverse lookup). + */ + virtual inline unsigned long GetEnthFromTemp(su2double& enthalpy, + const su2double val_prog, + const su2double val_temp, + su2double initial_value) { return 0; } + /*! * \brief Get fluid dynamic viscosity. */ @@ -319,7 +338,7 @@ class CFluidModel { * \brief Virtual member. * \param[in] T - Temperature value at the point. */ - virtual void SetTDState_T(su2double val_Temperature, const su2double* val_scalars = nullptr) {} + virtual void SetTDState_T(su2double val_Temperature, const su2double* val_scalars = nullptr) { } /*! * \brief Set fluid eddy viscosity provided by a turbulence model needed for computing effective thermal conductivity. diff --git a/SU2_CFD/include/numerics/species/species_diffusion.hpp b/SU2_CFD/include/numerics/species/species_diffusion.hpp index e100496b135..7db4d407cd7 100644 --- a/SU2_CFD/include/numerics/species/species_diffusion.hpp +++ b/SU2_CFD/include/numerics/species/species_diffusion.hpp @@ -71,8 +71,8 @@ class CAvgGrad_Species final : public CAvgGrad_Scalar { void FinishResidualCalc(const CConfig* config) override { for (auto iVar = 0u; iVar < nVar; iVar++) { - /* --- in case of species transport, Diffusion_Coeff is the binary diffusion coefficient --- */ const su2double Diffusivity_Lam = 0.5 * (Density_i * Diffusion_Coeff_i[iVar] + Density_j * Diffusion_Coeff_j[iVar]); + su2double Diffusivity_Turb = 0.0; if (turbulence) { diff --git a/SU2_CFD/include/numerics/species/species_sources.hpp b/SU2_CFD/include/numerics/species/species_sources.hpp index 7dd564720ef..87910f37d37 100644 --- a/SU2_CFD/include/numerics/species/species_sources.hpp +++ b/SU2_CFD/include/numerics/species/species_sources.hpp @@ -89,4 +89,5 @@ class CSourceAxisymmetric_Species : public CSourceBase_Species { * \return Lightweight const-view of residual and Jacobian. */ ResidualType<> ComputeResidual(const CConfig* config) override; + }; diff --git a/SU2_CFD/include/output/CFlowIncOutput.hpp b/SU2_CFD/include/output/CFlowIncOutput.hpp index bdb7a23c5cb..7a8621f25db 100644 --- a/SU2_CFD/include/output/CFlowIncOutput.hpp +++ b/SU2_CFD/include/output/CFlowIncOutput.hpp @@ -41,6 +41,7 @@ class CFlowIncOutput final: public CFlowOutput { TURB_MODEL turb_model; /*!< \brief The kind of turbulence model*/ bool heat; /*!< \brief Boolean indicating whether have a heat problem*/ bool weakly_coupled_heat; /*!< \brief Boolean indicating whether have a weakly coupled heat equation*/ + bool flamelet; /*!< \brief Boolean indicating whether we solve the flamelet equations */ unsigned short streamwisePeriodic; /*!< \brief Boolean indicating whether it is a streamwise periodic simulation. */ bool streamwisePeriodic_temperature; /*!< \brief Boolean indicating streamwise periodic temperature is used. */ diff --git a/SU2_CFD/include/output/CFlowOutput.hpp b/SU2_CFD/include/output/CFlowOutput.hpp index dcab274bf02..a47cae8de2b 100644 --- a/SU2_CFD/include/output/CFlowOutput.hpp +++ b/SU2_CFD/include/output/CFlowOutput.hpp @@ -125,7 +125,31 @@ class CFlowOutput : public CFVMOutput{ void SetVolumeOutputFieldsScalarResidual(const CConfig* config); /*! - * \brief Add scalar (turbulence/species) volume limiter fields (and more) for a point (FVMComp, FVMInc, FVMNEMO). + * \brief Add scalar (turbulence/species) volume primitive fields for a point (FVMComp, FVMInc, FVMNEMO). + * \param[in] config - Definition of the particular problem. + */ + void SetVolumeOutputFieldsScalarPrimitive(const CConfig* config); + + /*! + * \brief Add scalar (turbulence/species) volume limiter fields for a point (FVMComp, FVMInc, FVMNEMO). + * \param[in] config - Definition of the particular problem. + */ + void SetVolumeOutputFieldsScalarLimiter(const CConfig* config); + + /*! + * \brief Add flamelet volume source term fields for a point (FVMComp, FVMInc, FVMNEMO). + * \param[in] config - Definition of the particular problem. + */ + void SetVolumeOutputFieldsScalarSource(const CConfig* config); + + /*! + * \brief Add flamelet volume lookup value fields for a point (FVMComp, FVMInc, FVMNEMO). + * \param[in] config - Definition of the particular problem. + */ + void SetVolumeOutputFieldsScalarLookup(const CConfig* config); + + /*! + * \brief Add miscellaneous scalar volume fields for a point (FVMComp, FVMInc, FVMNEMO). * \param[in] config - Definition of the particular problem. */ void SetVolumeOutputFieldsScalarMisc(const CConfig* config); diff --git a/SU2_CFD/include/solvers/CSpeciesFlameletSolver.hpp b/SU2_CFD/include/solvers/CSpeciesFlameletSolver.hpp new file mode 100644 index 00000000000..95b4afe7311 --- /dev/null +++ b/SU2_CFD/include/solvers/CSpeciesFlameletSolver.hpp @@ -0,0 +1,171 @@ +/*! + * \file CSpeciesFlameletSolver.hpp + * \brief Headers of the CSpeciesFlameletSolver class + * \author D. Mayer, N. Beishuizen, T. Economon + * \version 7.5.1 "Blackbird" + * + * SU2 Project Website: https://su2code.github.io + * + * The SU2 Project is maintained by the SU2 Foundation + * (http://su2foundation.org) + * + * Copyright 2012-2023, SU2 Contributors (cf. AUTHORS.md) + * + * SU2 is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public + * License as published by the Free Software Foundation; either + * version 2.1 of the License, or (at your option) any later version. + * + * SU2 is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with SU2. If not, see . + */ + +#pragma once + +#include "CSpeciesSolver.hpp" + +/*! + * \class CSpeciesFlameletSolver + * \brief Main class for defining the flamelet model solver. + * \author N. Beishuizen + * \ingroup Scalar_Transport + */ +class CSpeciesFlameletSolver final : public CSpeciesSolver { + private: + vector conjugate_var; /*!< \brief CHT variables for each boundary and vertex. */ + + /*! + * \brief Compute the preconditioner for low-Mach flows. + * \param[in] geometry - Geometrical definition of the problem. + * \param[in] solver_container - Container vector with all the solutions. + * \param[in] config - Definition of the particular problem. + */ + void SetPreconditioner(CGeometry* geometry, CSolver** solver_container, CConfig* config); + + /*! + * \brief Compute the primitive variables (diffusivities). + * \param[in] solver_container - Container vector with all the solutions. + * \param[in] config - Definition of the particular problem. + * \param[in] Output - Boolean to determine whether to print output. + * \return - The number of non-physical points. + */ + unsigned long SetPrimitive_Variables(CSolver** solver_container, CConfig* config, bool Output); + + /*! + * \brief Generic implementation of the isothermal wall also covering CHT cases, + * for which the wall temperature is given by GetConjugateHeatVariable. + */ + void BC_Isothermal_Wall_Generic(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, + CNumerics* visc_numerics, CConfig* config, unsigned short val_marker, + bool cht_mode = false); + + public: + /*! + * \brief Constructor. + * \param[in] geometry - Geometrical definition of the problem. + * \param[in] config - Definition of the particular problem. + * \param[in] iMesh - Index of the mesh in multigrid computations. + * \param[in] FluidModel + */ + CSpeciesFlameletSolver(CGeometry* geometry, CConfig* config, unsigned short iMesh); + + /*! + * \brief Restart residual and compute gradients. + * \param[in] geometry - Geometrical definition of the problem. + * \param[in] solver_container - Container vector with all the solutions. + * \param[in] config - Definition of the particular problem. + * \param[in] iMesh - Index of the mesh in multigrid computations. + * \param[in] iRKStep - Current step of the Runge-Kutta iteration. + * \param[in] RunTime_EqSystem - System of equations which is going to be solved. + * \param[in] Output - Boolean to determine whether to print output. + */ + void Preprocessing(CGeometry* geometry, CSolver** solver_container, CConfig* config, unsigned short iMesh, + unsigned short iRKStep, unsigned short RunTime_EqSystem, bool Output) override; + + /*! + * \brief Set the initial condition for the scalar transport problem. + * \param[in] geometry - Geometrical definition of the problem. + * \param[in] solver_container - Container with all the solutions. + * \param[in] config - Definition of the particular problem. + * \param[in] ExtIter - External iteration. + */ + void SetInitialCondition(CGeometry** geometry, CSolver*** solver_container, CConfig* config, + unsigned long ExtIter) override; + + /*! + * \brief Source term computation. + * \param[in] geometry - Geometrical definition of the problem. + * \param[in] solver_container - Container vector with all the solutions. + * \param[in] numerics - Description of the numerical method. + * \param[in] config - Definition of the particular problem. + * \param[in] iMesh - Index of the mesh in multigrid computations. + */ + void Source_Residual(CGeometry* geometry, CSolver** solver_container, CNumerics** numerics_container, CConfig* config, + unsigned short iMesh) override; + + /*! + * \brief Impose the Navier-Stokes 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 inlet 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_Inlet(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, CNumerics* visc_numerics, + CConfig* config, unsigned short val_marker) override; + + /*! + * \brief Impose the (received) conjugate heat variables. + * \param[in] geometry - Geometrical definition of the problem. + * \param[in] solver_container - Container vector with all the solutions. + * \param[in] 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_ConjugateHeat_Interface(CGeometry* geometry, CSolver** solver_container, CNumerics* numerics, CConfig* config, + unsigned short val_marker) override; + + /*! + * \brief Get the conjugate heat variables. + * \param[in] val_marker - The marker index. + * \param[in] val_vertex - The vertex index. + * \param[in] pos_var - The variable position (in vector of all conjugate heat variables). + */ + inline su2double GetConjugateHeatVariable(unsigned short val_marker, unsigned long val_vertex, + unsigned short pos_var) const override { + return conjugate_var[val_marker][val_vertex][pos_var]; + } + + /*! + * \brief Set the conjugate heat variables. + * \param[in] val_marker - The marker index. + * \param[in] val_vertex - The vertex index. + * \param[in] pos_var - The variable position (in vector of all conjugate heat variables). + * \param[in] relaxation_factor - The relaxation factor for the change of the variables. + * \param[in] val_var - The value of the variable. + */ + inline void SetConjugateHeatVariable(unsigned short val_marker, unsigned long val_vertex, unsigned short pos_var, + su2double relaxation_factor, su2double val_var) override { + conjugate_var[val_marker][val_vertex][pos_var] = + relaxation_factor * val_var + (1.0 - relaxation_factor) * conjugate_var[val_marker][val_vertex][pos_var]; + } + +}; diff --git a/SU2_CFD/include/solvers/CSpeciesSolver.hpp b/SU2_CFD/include/solvers/CSpeciesSolver.hpp index d4d751e3ff5..1d4eef366de 100644 --- a/SU2_CFD/include/solvers/CSpeciesSolver.hpp +++ b/SU2_CFD/include/solvers/CSpeciesSolver.hpp @@ -65,6 +65,8 @@ class CSpeciesSolver : public CScalarSolver { */ void LoadRestart(CGeometry** geometry, CSolver*** solver, CConfig* config, int val_iter, bool val_update_geo) final; + void Initialize(CGeometry* geometry, CConfig* config, unsigned short iMesh, unsigned short nVar); + /*! * \brief Restart residual and compute gradients. * \param[in] geometry - Geometrical definition of the problem. @@ -76,7 +78,7 @@ class CSpeciesSolver : public CScalarSolver { * \param[in] Output - boolean to determine whether to print output. */ void Preprocessing(CGeometry* geometry, CSolver** solver_container, CConfig* config, unsigned short iMesh, - unsigned short iRKStep, unsigned short RunTime_EqSystem, bool Output) final; + unsigned short iRKStep, unsigned short RunTime_EqSystem, bool Output) override; /*! * \brief Compute the viscous flux for the turbulent equation at a particular edge. @@ -100,7 +102,7 @@ class CSpeciesSolver : public CScalarSolver { * \param[in] val_marker - Surface marker where the boundary condition is applied. */ void BC_Inlet(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, CNumerics* visc_numerics, - CConfig* config, unsigned short val_marker) final; + CConfig* config, unsigned short val_marker) override; /*! * \brief Store of a set of provided inlet profile values at a vertex. diff --git a/SU2_CFD/include/variables/CSpeciesFlameletVariable.hpp b/SU2_CFD/include/variables/CSpeciesFlameletVariable.hpp new file mode 100644 index 00000000000..6f2eb6fe42b --- /dev/null +++ b/SU2_CFD/include/variables/CSpeciesFlameletVariable.hpp @@ -0,0 +1,90 @@ +/*! + * \file CSpeciesFlameletVariable.hpp + * \brief Base class for defining the variables of the flamelet transport model. + * \author D. Mayer, T. Economon, N. Beishuizen + * \version 7.5.1 "Blackbird" + * + * SU2 Project Website: https://su2code.github.io + * + * The SU2 Project is maintained by the SU2 Foundation + * (http://su2foundation.org) + * + * Copyright 2012-2023, SU2 Contributors (cf. AUTHORS.md) + * + * SU2 is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public + * License as published by the Free Software Foundation; either + * version 2.1 of the License, or (at your option) any later version. + * + * SU2 is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with SU2. If not, see . + */ + +#pragma once + +#include "CSpeciesVariable.hpp" + +/*! + * \class CSpeciesFlameletVariable + * \brief Base class for defining the variables of the flamelet model. + */ +class CSpeciesFlameletVariable final : public CSpeciesVariable { + protected: + MatrixType source_scalar; /*!< \brief Vector of the source terms from the lookup table for each scalar equation */ + MatrixType lookup_scalar; /*!< \brief Vector of the source terms from the lookup table for each scalar equation */ + su2vector table_misses; /*!< \brief Vector of lookup table misses. */ + + public: + /*! + * \brief Constructor of the class. + * \param[in] species_inf - species variable values (initialization value). + * \param[in] npoint - Number of points/nodes/vertices in the domain. + * \param[in] ndim - Number of dimensions of the problem. + * \param[in] nvar - Number of variables of the problem. + * \param[in] config - Definition of the particular problem. + */ + CSpeciesFlameletVariable(const su2double* species_inf, unsigned long npoint, unsigned long ndim, unsigned long nvar, + const CConfig* config); + + /*! + * \brief Set the value of the transported scalar source term. + * \param[in] iPoint - the location where the value has to be set. + * \param[in] val_lookup_scalar - the value of the scalar to set. + * \param[in] val_ivar - Eqn. index to the transport equation. + */ + inline void SetLookupScalar(unsigned long iPoint, su2double val_lookup_scalar, unsigned short val_ivar) override { + lookup_scalar(iPoint, val_ivar) = val_lookup_scalar; + } + + /*! + + * \brief Set a source term for the specie transport equation. + * \param[in] iPoint - Node index. + * \param[in] val_ivar - Species index. + * \param[in] val_source - Source term value. + */ + inline void SetScalarSource(unsigned long iPoint, unsigned short val_ivar, su2double val_source) override { + source_scalar(iPoint, val_ivar) = val_source; + } + + /*! + * \brief Get the value of the transported scalars source term. + * \return Pointer to the transported scalars source term. + */ + inline const su2double* GetScalarSources(unsigned long iPoint) const override { return source_scalar[iPoint]; } + + /*! + * \brief Get the value of the looked up table based on the transported scalar. + * \return Pointer to the transported scalars source term. + */ + inline const su2double* GetScalarLookups(unsigned long iPoint) const override { return lookup_scalar[iPoint]; } + + inline void SetTableMisses(unsigned long iPoint, unsigned short misses) override { table_misses[iPoint] = misses; } + + inline unsigned short GetTableMisses(unsigned long iPoint) const override { return table_misses[iPoint]; } +}; diff --git a/SU2_CFD/include/variables/CSpeciesVariable.hpp b/SU2_CFD/include/variables/CSpeciesVariable.hpp index d668e006db2..a3d254ec30b 100644 --- a/SU2_CFD/include/variables/CSpeciesVariable.hpp +++ b/SU2_CFD/include/variables/CSpeciesVariable.hpp @@ -38,7 +38,7 @@ class CSpeciesVariable : public CScalarVariable { MatrixType Diffusivity; /*!< \brief Matrix (nPoint,nVar) of mass diffusivities for scalar transport. */ public: - static constexpr size_t MAXNVAR = 4; /*!< \brief Max number of variables for static arrays. Increase, if necessary. */ + static constexpr size_t MAXNVAR = 20; /*!< \brief Max number of variables for static arrays. Increase, if necessary. */ /*! * \brief Constructor of the class. diff --git a/SU2_CFD/include/variables/CVariable.hpp b/SU2_CFD/include/variables/CVariable.hpp index a5beafd75e0..5167ef953f6 100644 --- a/SU2_CFD/include/variables/CVariable.hpp +++ b/SU2_CFD/include/variables/CVariable.hpp @@ -2312,4 +2312,17 @@ class CVariable { virtual su2double GetSourceTerm_DispAdjoint(unsigned long iPoint, unsigned long iDim) const { return 0.0; } virtual su2double GetSourceTerm_VelAdjoint(unsigned long iPoint, unsigned long iDim) const { return 0.0; } + /*! + * \brief LUT premixed flamelet: virtual functions for the speciesflameletvariable LUT + */ + inline virtual void SetLookupScalar(unsigned long iPoint, su2double val_lookup_scalar, unsigned short val_ivar) { } + + inline virtual void SetScalarSource(unsigned long iPoint, unsigned short val_ivar, su2double val_source) { } + + inline virtual void SetTableMisses(unsigned long iPoint, unsigned short misses) { } + + inline virtual unsigned short GetTableMisses(unsigned long iPoint) const { return 0; } + + inline virtual const su2double *GetScalarSources(unsigned long iPoint) const { return nullptr; } + inline virtual const su2double *GetScalarLookups(unsigned long iPoint) const { return nullptr; } }; diff --git a/SU2_CFD/src/drivers/CDriver.cpp b/SU2_CFD/src/drivers/CDriver.cpp index 639e34c4304..1eee1af5e85 100644 --- a/SU2_CFD/src/drivers/CDriver.cpp +++ b/SU2_CFD/src/drivers/CDriver.cpp @@ -1430,6 +1430,7 @@ void CDriver::InstantiateSpeciesNumerics(unsigned short nVar_Species, int offset numerics[iMGlevel][SPECIES_SOL][source_second_term] = new CSourceNothing(nDim, nVar_Species, config); } } + /*--- Explicit instantiation of the template above, needed because it is defined in a cpp file, instead of hpp. ---*/ template void CDriver::InstantiateSpeciesNumerics>( unsigned short, int, const CConfig*, const CSolver*, CNumerics****&) const; @@ -3034,6 +3035,8 @@ void CFluidDriver::Preprocess(unsigned long Iter) { if (config_container[iZone]->GetFluidProblem()) { for (iInst = 0; iInst < nInst[iZone]; iInst++) { solver_container[iZone][iInst][MESH_0][FLOW_SOL]->SetInitialCondition(geometry_container[iZone][INST_0], solver_container[iZone][iInst], config_container[iZone], Iter); + if (config_container[iZone]->GetKind_Species_Model() != SPECIES_MODEL::NONE) + solver_container[iZone][iInst][MESH_0][SPECIES_SOL]->SetInitialCondition(geometry_container[iZone][INST_0], solver_container[iZone][iInst], config_container[iZone], Iter); } } } diff --git a/SU2_CFD/src/drivers/CMultizoneDriver.cpp b/SU2_CFD/src/drivers/CMultizoneDriver.cpp index 61f59223fcf..b3a7af5e369 100644 --- a/SU2_CFD/src/drivers/CMultizoneDriver.cpp +++ b/SU2_CFD/src/drivers/CMultizoneDriver.cpp @@ -245,6 +245,14 @@ void CMultizoneDriver::Preprocess(unsigned long TimeIter) { solver_container[iZone][INST_0], config_container[iZone], TimeIter); } + if (!fsi && (config_container[iZone]->GetKind_Species_Model() != SPECIES_MODEL::NONE)) { + /*--- Set the initial condition for SPECIES solver (species or flamelet) ----------------------*/ + solver_container[iZone][INST_0][MESH_0][SPECIES_SOL]->SetInitialCondition(geometry_container[ZONE_0][INST_0], + solver_container[ZONE_0][INST_0], + config_container[ZONE_0], TimeIter); + } + + } /*--- Ramp turbo values for unsteady problems here, otherwise do it over outer iterations. ---*/ diff --git a/SU2_CFD/src/drivers/CSinglezoneDriver.cpp b/SU2_CFD/src/drivers/CSinglezoneDriver.cpp index 964888a83fb..c1f2c2519c8 100644 --- a/SU2_CFD/src/drivers/CSinglezoneDriver.cpp +++ b/SU2_CFD/src/drivers/CSinglezoneDriver.cpp @@ -134,6 +134,11 @@ void CSinglezoneDriver::Preprocess(unsigned long TimeIter) { solver_container[ZONE_0][INST_0], config_container[ZONE_0], TimeIter); } + if (config_container[ZONE_0]->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + solver_container[ZONE_0][INST_0][MESH_0][SPECIES_SOL]->SetInitialCondition(geometry_container[ZONE_0][INST_0], + solver_container[ZONE_0][INST_0], + config_container[ZONE_0], TimeIter); + } else if (config_container[ZONE_0]->GetHeatProblem()) { /*--- Set the initial condition for HEAT equation ---------------------------------------------*/ solver_container[ZONE_0][INST_0][MESH_0][HEAT_SOL]->SetInitialCondition(geometry_container[ZONE_0][INST_0], diff --git a/SU2_CFD/src/fluid/CFluidFlamelet.cpp b/SU2_CFD/src/fluid/CFluidFlamelet.cpp new file mode 100644 index 00000000000..34051288559 --- /dev/null +++ b/SU2_CFD/src/fluid/CFluidFlamelet.cpp @@ -0,0 +1,138 @@ +/*! + * \file CfluidFlamelet.cpp + * \brief Main subroutines of CFluidFlamelet class + * \author D. Mayer, T. Economon, N. Beishuizen + * \version 7.5.1 "Blackbird" + * + * SU2 Project Website: https://su2code.github.io + * + * The SU2 Project is maintained by the SU2 Foundation + * (http://su2foundation.org) + * + * Copyright 2012-2023, SU2 Contributors (cf. AUTHORS.md) + * + * SU2 is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public + * License as published by the Free Software Foundation; either + * version 2.1 of the License, or (at your option) any later version. + * + * SU2 is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with SU2. If not, see . + */ + +#include "../include/fluid/CFluidFlamelet.hpp" +#include "../../../Common/include/containers/CLookUpTable.hpp" + +CFluidFlamelet::CFluidFlamelet(CConfig* config, su2double value_pressure_operating) : CFluidModel() { + rank = SU2_MPI::GetRank(); + + /* -- number of auxiliary species transport equations, e.g. 1=CO, 2=NOx --- */ + n_user_scalars = config->GetNUserScalars(); + n_control_vars = config->GetNControlVars(); + n_scalars = config->GetNScalars(); + + if (rank == MASTER_NODE) { + cout << "Number of scalars: " << n_scalars << endl; + cout << "Number of user scalars: " << n_user_scalars << endl; + cout << "Number of control variables: " << n_control_vars << endl; + } + + if (rank == MASTER_NODE) { + cout << "*****************************************" << endl; + cout << "*** initializing the lookup table ***" << endl; + cout << "*****************************************" << endl; + } + + table_scalar_names.resize(n_scalars); + table_scalar_names[I_ENTH] = "EnthalpyTot"; + table_scalar_names[I_PROGVAR] = "ProgressVariable"; + /*--- auxiliary species transport equations---*/ + for (size_t i_aux = 0; i_aux < n_user_scalars; i_aux++) { + table_scalar_names[n_control_vars + i_aux] = config->GetUserScalarName(i_aux); + } + + look_up_table = new CLookUpTable(config->GetFileNameLUT(), table_scalar_names[I_PROGVAR], table_scalar_names[I_ENTH]); + + n_lookups = config->GetNLookups(); + table_lookup_names.resize(n_lookups); + for (int i_lookup = 0; i_lookup < n_lookups; ++i_lookup) { + table_lookup_names[i_lookup] = config->GetLUTLookupName(i_lookup); + } + + Pressure = value_pressure_operating; + + /*--- Thermodynamic state variables and names. ---*/ + varnames_TD.resize(LOOKUP_TD::SIZE); + val_vars_TD.resize(LOOKUP_TD::SIZE); + + /*--- The string in varnames_TD as it appears in the LUT file. ---*/ + varnames_TD[LOOKUP_TD::TEMPERATURE] = "Temperature"; + varnames_TD[LOOKUP_TD::DENSITY] = "Density"; + varnames_TD[LOOKUP_TD::HEATCAPACITY] = "Cp"; + varnames_TD[LOOKUP_TD::VISCOSITY] = "ViscosityDyn"; + varnames_TD[LOOKUP_TD::CONDUCTIVITY] = "Conductivity"; + varnames_TD[LOOKUP_TD::DIFFUSIONCOEFFICIENT] = "DiffusionCoefficient"; + varnames_TD[LOOKUP_TD::MOLARWEIGHT] = "MolarWeightMix"; + +} + +CFluidFlamelet::~CFluidFlamelet() { delete look_up_table; } + +void CFluidFlamelet::SetTDState_T(su2double val_temperature, const su2double* val_scalars) { + su2double val_enth = val_scalars[I_ENTH]; + su2double val_prog = val_scalars[I_PROGVAR]; + + /*--- Add all quantities and their names to the look up vectors. ---*/ + look_up_table->LookUp_XY(varnames_TD, val_vars_TD, val_prog, val_enth); + + Temperature = val_vars_TD[LOOKUP_TD::TEMPERATURE]; + Density = val_vars_TD[LOOKUP_TD::DENSITY]; + Cp = val_vars_TD[LOOKUP_TD::HEATCAPACITY]; + Mu = val_vars_TD[LOOKUP_TD::VISCOSITY]; + Kt = val_vars_TD[LOOKUP_TD::CONDUCTIVITY]; + mass_diffusivity = val_vars_TD[LOOKUP_TD::DIFFUSIONCOEFFICIENT]; + molar_weight = val_vars_TD[LOOKUP_TD::MOLARWEIGHT]; + + /*--- Compute Cv from Cp and molar weight of the mixture (ideal gas). ---*/ + Cv = Cp - UNIVERSAL_GAS_CONSTANT / molar_weight; +} + +/* --- Total enthalpy is the transported variable, but we usually have temperature as a boundary condition, + so we do a reverse lookup */ +unsigned long CFluidFlamelet::GetEnthFromTemp(su2double& val_enth, const su2double val_prog, const su2double val_temp, + const su2double initial_value) { + /*--- convergence criterion for temperature in [K], high accuracy needed for restarts. ---*/ + su2double delta_temp_final = 0.001; + su2double enth_iter = initial_value; + su2double delta_enth; + su2double delta_temp_iter = 1e10; + unsigned long exit_code = 0; + const int counter_limit = 1000; + + int counter = 0; + while ((abs(delta_temp_iter) > delta_temp_final) && (counter++ < counter_limit)) { + /*--- Add all quantities and their names to the look up vectors. ---*/ + look_up_table->LookUp_XY(varnames_TD, val_vars_TD, val_prog, enth_iter); + Temperature = val_vars_TD[LOOKUP_TD::TEMPERATURE]; + Cp = val_vars_TD[LOOKUP_TD::HEATCAPACITY]; + + delta_temp_iter = val_temp - Temperature; + + delta_enth = Cp * delta_temp_iter; + + enth_iter += delta_enth; + } + + val_enth = enth_iter; + + if (counter >= counter_limit) { + exit_code = 1; + } + + return exit_code; +} diff --git a/SU2_CFD/src/fluid/CFluidModel.cpp b/SU2_CFD/src/fluid/CFluidModel.cpp index 15ab39b24e2..5aca9e783da 100644 --- a/SU2_CFD/src/fluid/CFluidModel.cpp +++ b/SU2_CFD/src/fluid/CFluidModel.cpp @@ -41,6 +41,7 @@ #include "../../include/fluid/CPolynomialConductivityRANS.hpp" #include "../../include/fluid/CPolynomialViscosity.hpp" #include "../../include/fluid/CSutherland.hpp" +#include "../../include/fluid/CFluidFlamelet.hpp" #include "../../include/fluid/CCoolPropViscosity.hpp" #include "../../include/fluid/CConstantLewisDiffusivity.hpp" #include "../../include/fluid/CCoolPropConductivity.hpp" @@ -58,6 +59,9 @@ unique_ptr CFluidModel::MakeLaminarViscosityModel(const CConfig case VISCOSITYMODEL::POLYNOMIAL: return unique_ptr>( new CPolynomialViscosity(config->GetMu_PolyCoeffND())); + case VISCOSITYMODEL::FLAMELET: + /*--- Viscosity is obtained from the LUT ---*/ + return nullptr; default: SU2_MPI::Error("Viscosity model not available.", CURRENT_FUNCTION); return nullptr; @@ -108,10 +112,12 @@ unique_ptr CFluidModel::MakeThermalConductivityModel(const C new CPolynomialConductivity(config->GetKt_PolyCoeffND())); } break; + case CONDUCTIVITYMODEL::FLAMELET: + /*--- Conductivity is obtained from the LUT ---*/ + return nullptr; default: SU2_MPI::Error("Conductivity model not available.", CURRENT_FUNCTION); return nullptr; - break; } } @@ -134,6 +140,10 @@ unique_ptr CFluidModel::MakeMassDiffusivityModel(const CConfi return unique_ptr( new CConstantLewisDiffusivity(config->GetConstant_Lewis_Number(iSpecies))); break; + case DIFFUSIVITYMODEL::FLAMELET: + /*--- Diffusivity is obtained from the LUT ---*/ + return nullptr; + break; default: SU2_MPI::Error("Diffusivity model not available.", CURRENT_FUNCTION); return nullptr; diff --git a/SU2_CFD/src/integration/CIntegration.cpp b/SU2_CFD/src/integration/CIntegration.cpp index 1da64f36a32..d27be847c53 100644 --- a/SU2_CFD/src/integration/CIntegration.cpp +++ b/SU2_CFD/src/integration/CIntegration.cpp @@ -161,7 +161,12 @@ void CIntegration::Space_Integration(CGeometry *geometry, solver_container[MainSolver]->BC_Custom(geometry, solver_container, conv_bound_numerics, visc_bound_numerics, config, iMarker); break; case CHT_WALL_INTERFACE: - if ((MainSolver == HEAT_SOL) || ((MainSolver == FLOW_SOL) && ((config->GetKind_Regime() == ENUM_REGIME::COMPRESSIBLE) || config->GetEnergy_Equation()))) { + + if ( MainSolver == FLOW_SOL && (config->GetKind_FluidModel() == FLUID_FLAMELET) ){ + solver_container[MainSolver]->BC_Isothermal_Wall(geometry, solver_container, conv_bound_numerics, visc_bound_numerics, config, iMarker); + break; + } + if ((MainSolver == SPECIES_SOL) || (MainSolver == HEAT_SOL) || ((MainSolver == FLOW_SOL) && ((config->GetKind_Regime() == ENUM_REGIME::COMPRESSIBLE) || config->GetEnergy_Equation()))) { solver_container[MainSolver]->BC_ConjugateHeat_Interface(geometry, solver_container, conv_bound_numerics, config, iMarker); } else { diff --git a/SU2_CFD/src/interfaces/cht/CConjugateHeatInterface.cpp b/SU2_CFD/src/interfaces/cht/CConjugateHeatInterface.cpp index 56fcea263e2..594863bb9a1 100644 --- a/SU2_CFD/src/interfaces/cht/CConjugateHeatInterface.cpp +++ b/SU2_CFD/src/interfaces/cht/CConjugateHeatInterface.cpp @@ -104,6 +104,7 @@ void CConjugateHeatInterface::GetDonor_Variable(CSolver *donor_solution, CGeomet switch (donor_config->GetKind_ConductivityModel()) { case CONDUCTIVITYMODEL::CONSTANT: + case CONDUCTIVITYMODEL::FLAMELET: case CONDUCTIVITYMODEL::COOLPROP: thermal_conductivity = thermal_conductivityND*donor_config->GetThermal_Conductivity_Ref(); break; diff --git a/SU2_CFD/src/iteration/CFluidIteration.cpp b/SU2_CFD/src/iteration/CFluidIteration.cpp index 298ee40d56c..da63e7812a5 100644 --- a/SU2_CFD/src/iteration/CFluidIteration.cpp +++ b/SU2_CFD/src/iteration/CFluidIteration.cpp @@ -98,7 +98,7 @@ void CFluidIteration::Iterate(COutput* output, CIntegration**** integration, CGe RUNTIME_TURB_SYS, val_iZone, val_iInst); } - if (config[val_iZone]->GetKind_Species_Model() != SPECIES_MODEL::NONE){ + if (config[val_iZone]->GetKind_Species_Model() != SPECIES_MODEL::NONE) { config[val_iZone]->SetGlobalParam(main_solver, RUNTIME_SPECIES_SYS); integration[val_iZone][val_iInst][SPECIES_SOL]->SingleGrid_Iteration(geometry, solver, numerics, config, RUNTIME_SPECIES_SYS, val_iZone, val_iInst); diff --git a/SU2_CFD/src/meson.build b/SU2_CFD/src/meson.build index 9ca0188ffc0..98f71f9fd45 100644 --- a/SU2_CFD/src/meson.build +++ b/SU2_CFD/src/meson.build @@ -6,6 +6,7 @@ su2_cfd_src = files(['definition_structure.cpp', su2_cfd_src += files(['fluid/CFluidModel.cpp', 'fluid/CIdealGas.cpp', 'fluid/CFluidScalar.cpp', + 'fluid/CFluidFlamelet.cpp', 'fluid/CPengRobinson.cpp', 'fluid/CVanDerWaalsGas.cpp', 'fluid/CCoolProp.cpp', @@ -58,6 +59,7 @@ su2_cfd_src += files(['variables/CIncNSVariable.cpp', 'variables/CTurbVariable.cpp', 'variables/CScalarVariable.cpp', 'variables/CSpeciesVariable.cpp', + 'variables/CSpeciesFlameletVariable.cpp', 'variables/CAdjNSVariable.cpp', 'variables/CBaselineVariable.cpp', 'variables/CDiscAdjFEABoundVariable.cpp', @@ -107,6 +109,7 @@ su2_cfd_src += files(['solvers/CSolverFactory.cpp', 'solvers/CSolver.cpp', 'solvers/CTemplateSolver.cpp', 'solvers/CSpeciesSolver.cpp', + 'solvers/CSpeciesFlameletSolver.cpp', 'solvers/CTransLMSolver.cpp', 'solvers/CTurbSolver.cpp', 'solvers/CTurbSASolver.cpp', diff --git a/SU2_CFD/src/numerics/species/species_sources.cpp b/SU2_CFD/src/numerics/species/species_sources.cpp index 06a26902b8a..ab51d4b15e7 100644 --- a/SU2_CFD/src/numerics/species/species_sources.cpp +++ b/SU2_CFD/src/numerics/species/species_sources.cpp @@ -90,7 +90,6 @@ CNumerics::ResidualType<> CSourceAxisymmetric_Species::ComputeResidual(const /*--- Contribution due to 2D axisymmetric formulation ---*/ if (Coord_i[1] > EPS) { - AD::SetPreaccIn(Coord_i[1]); AD::SetPreaccIn(Diffusion_Coeff_i, nVar); AD::SetPreaccIn(ScalarVar_Grad_i, nVar, nDim); @@ -123,7 +122,7 @@ CNumerics::ResidualType<> CSourceAxisymmetric_Species::ComputeResidual(const if (turbulence) Mass_Diffusivity_Tur = V_i[idx.EddyViscosity()] / Sc_t; - for (auto iVar=0u; iVar < nVar; iVar++){ + for (auto iVar = 0u; iVar < nVar; iVar++) { residual[iVar] += yinv * Volume * (Density_i * Diffusion_Coeff_i[iVar] + Mass_Diffusivity_Tur) * ScalarVar_Grad_i[iVar][1]; } } diff --git a/SU2_CFD/src/output/CAdjFlowOutput.cpp b/SU2_CFD/src/output/CAdjFlowOutput.cpp index 4ef01b6c16c..9489926a537 100644 --- a/SU2_CFD/src/output/CAdjFlowOutput.cpp +++ b/SU2_CFD/src/output/CAdjFlowOutput.cpp @@ -55,11 +55,22 @@ void CAdjFlowOutput::AddHistoryOutputFields_AdjScalarRMS_RES(const CConfig* conf } } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) { for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { AddHistoryOutput("RMS_ADJ_SPECIES_" + std::to_string(iVar), "rms[A_rho*Y_" + std::to_string(iVar) + "]", ScreenOutputFormat::FIXED, "RMS_RES", "Root-mean square residual of the adjoint transported species.", HistoryFieldType::RESIDUAL); } } + + if (config->GetKind_Species_Model() == SPECIES_MODEL::FLAMELET) { + + AddHistoryOutput("RMS_ADJ_PROGRESS_VARIABLE", "rms[PV]", ScreenOutputFormat::FIXED, "RMS_RES", "Root-mean square residual of the adjoint progress variable.", HistoryFieldType::RESIDUAL); + AddHistoryOutput("RMS_ADJ_TOTAL_ENTHALPY", "rms[Enth]", ScreenOutputFormat::FIXED, "RMS_RES", "Root-mean square residual of the adjoint total enthalpy.", HistoryFieldType::RESIDUAL); + + for (unsigned short i_scalar = 0; i_scalar < config->GetNUserScalars(); i_scalar++) { + const auto& scalar_name = config->GetUserScalarName(i_scalar); + AddHistoryOutput("RMS_ADJ_" + scalar_name, "rms[" + scalar_name + "]", ScreenOutputFormat::FIXED, "RMS_RES", "Root-mean square residual of the adjoint of " + scalar_name + " .", HistoryFieldType::RESIDUAL); + } + } } void CAdjFlowOutput::AddHistoryOutputFields_AdjScalarMAX_RES(const CConfig* config) { @@ -80,11 +91,22 @@ void CAdjFlowOutput::AddHistoryOutputFields_AdjScalarMAX_RES(const CConfig* conf } } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) { for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { AddHistoryOutput("MAX_ADJ_SPECIES_" + std::to_string(iVar), "max[A_rho*Y_" + std::to_string(iVar) + "]",ScreenOutputFormat::FIXED, "MAX_RES", "Maximum residual of the adjoint transported species.", HistoryFieldType::RESIDUAL); } } + + if (config->GetKind_Species_Model() == SPECIES_MODEL::FLAMELET) { + + AddHistoryOutput("MAX_ADJ_PROGRESS_VARIABLE", "max[PV]",ScreenOutputFormat::FIXED, "MAX_RES", "Maximum residual of the adjoint progress variable.", HistoryFieldType::RESIDUAL); + AddHistoryOutput("MAX_ADJ_TOTAL_ENTHALPY", "max[Enth]",ScreenOutputFormat::FIXED, "MAX_RES", "Maximum residual of the adjoint total enthalpy.", HistoryFieldType::RESIDUAL); + + for (unsigned short i_scalar = 0; i_scalar < config->GetNUserScalars(); i_scalar++) { + const auto& scalar_name = config->GetUserScalarName(i_scalar); + AddHistoryOutput("MAX_ADJ_" + scalar_name, "max[scalar_" + scalar_name + "]",ScreenOutputFormat::FIXED, "MAX_RES", "Maximum residual of the adjoint of " + scalar_name + " .", HistoryFieldType::RESIDUAL); + } + } } void CAdjFlowOutput::AddHistoryOutputFields_AdjScalarBGS_RES(const CConfig* config) { @@ -107,11 +129,22 @@ void CAdjFlowOutput::AddHistoryOutputFields_AdjScalarBGS_RES(const CConfig* conf } } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) { for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { AddHistoryOutput("BGS_ADJ_SPECIES_" + std::to_string(iVar), "bgs[A_rho*Y_" + std::to_string(iVar) + "]", ScreenOutputFormat::FIXED, "BGS_RES", "BGS residual of the adjoint transported species.", HistoryFieldType::RESIDUAL); } } + + if (config->GetKind_Species_Model() == SPECIES_MODEL::FLAMELET) { + + AddHistoryOutput("BGS_ADJ_PROGRESS_VARIABLE", "bgs[PV]", ScreenOutputFormat::FIXED, "BGS_RES", "BGS residual of the adjoint progress variable.", HistoryFieldType::RESIDUAL); + AddHistoryOutput("BGS_ADJ_TOTAL_ENTHALPY", "bgs[Enth]", ScreenOutputFormat::FIXED, "BGS_RES", "BGS residual of the adjoint total enthalpy.", HistoryFieldType::RESIDUAL); + + for (unsigned short i_scalar = 0; i_scalar < config->GetNUserScalars(); i_scalar++) { + const auto& scalar_name = config->GetUserScalarName(i_scalar); + AddHistoryOutput("BGS_ADJ_" + scalar_name, "bgs[" + scalar_name + "]", ScreenOutputFormat::FIXED, "BGS_RES", "BGS residual of the adjoint of " + scalar_name + " .", HistoryFieldType::RESIDUAL); + } + } } void CAdjFlowOutput::AddHistoryOutputFieldsAdjScalarLinsol(const CConfig* config) { @@ -120,10 +153,15 @@ void CAdjFlowOutput::AddHistoryOutputFieldsAdjScalarLinsol(const CConfig* config AddHistoryOutput("LINSOL_RESIDUAL_TURB", "LinSolResTurb", ScreenOutputFormat::FIXED, "LINSOL", "Residual of the linear solver for turbulence."); } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) { AddHistoryOutput("LINSOL_ITER_SPECIES", "LinSolIterSpecies", ScreenOutputFormat::INTEGER, "LINSOL", "Number of iterations of the linear solver for species solver."); AddHistoryOutput("LINSOL_RESIDUAL_SPECIES", "LinSolResSpecies", ScreenOutputFormat::FIXED, "LINSOL", "Residual of the linear solver for species solver."); } + + if (config->GetKind_Species_Model() == SPECIES_MODEL::FLAMELET) { + AddHistoryOutput("LINSOL_ITER_FLAMELET", "LinSolIterScalar", ScreenOutputFormat::INTEGER, "LINSOL", "Number of iterations of the linear solver for scalar solver."); + AddHistoryOutput("LINSOL_RESIDUAL_FLAMELET", "LinSolResScalar", ScreenOutputFormat::FIXED, "LINSOL", "Residual of the linear solver for scalar solver."); + } } // clang-format on @@ -160,7 +198,7 @@ void CAdjFlowOutput::LoadHistoryDataAdjScalar(const CConfig* config, const CSolv } } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) { for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { SetHistoryOutputValue("RMS_ADJ_SPECIES_" + std::to_string(iVar), log10(adjspecies_solver->GetRes_RMS(iVar))); SetHistoryOutputValue("MAX_ADJ_SPECIES_" + std::to_string(iVar), log10(adjspecies_solver->GetRes_Max(iVar))); @@ -173,6 +211,30 @@ void CAdjFlowOutput::LoadHistoryDataAdjScalar(const CConfig* config, const CSolv SetHistoryOutputValue("LINSOL_RESIDUAL_SPECIES", log10(adjspecies_solver->GetResLinSolver())); } + if (config->GetKind_Species_Model() == SPECIES_MODEL::FLAMELET) { + + SetHistoryOutputValue("RMS_ADJ_PROGRESS_VARIABLE", log10(adjspecies_solver->GetRes_RMS(I_PROGVAR))); + SetHistoryOutputValue("MAX_ADJ_PROGRESS_VARIABLE", log10(adjspecies_solver->GetRes_Max(I_PROGVAR))); + SetHistoryOutputValue("RMS_ADJ_TOTAL_ENTHALPY", log10(adjspecies_solver->GetRes_RMS(I_ENTH))); + SetHistoryOutputValue("MAX_ADJ_TOTAL_ENTHALPY", log10(adjspecies_solver->GetRes_Max(I_ENTH))); + if (multiZone) { + SetHistoryOutputValue("BGS_ADJ_PROGRESS_VARIABLE", log10(adjspecies_solver->GetRes_BGS(I_PROGVAR))); + SetHistoryOutputValue("BGS_ADJ_TOTAL_ENTHALPY", log10(adjspecies_solver->GetRes_BGS(I_ENTH))); + } + + for (unsigned short i_scalar = 0; i_scalar < config->GetNUserScalars(); i_scalar++) { + const auto& scalar_name = config->GetUserScalarName(i_scalar); + SetHistoryOutputValue("RMS_ADJ_" + scalar_name, log10(adjspecies_solver->GetRes_RMS(2 + i_scalar))); + SetHistoryOutputValue("MAX_ADJ_" + scalar_name, log10(adjspecies_solver->GetRes_Max(2 + i_scalar))); + if (multiZone) { + SetHistoryOutputValue("BGS_ADJ_" + scalar_name, log10(adjspecies_solver->GetRes_BGS(2 + i_scalar))); + } + } + + SetHistoryOutputValue("LINSOL_ITER_FLAMELET", adjspecies_solver->GetIterLinSolver()); + SetHistoryOutputValue("LINSOL_RESIDUAL_FLAMELET", log10(adjspecies_solver->GetResLinSolver())); + } + ComputeSimpleCustomOutputs(config); } @@ -194,12 +256,23 @@ void CAdjFlowOutput::SetVolumeOutputFieldsAdjScalarSolution(const CConfig* confi } } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) { for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { AddVolumeOutput("ADJ_SPECIES_" + std::to_string(iVar), "Adjoint_Species_" + std::to_string(iVar), "SOLUTION", "Adjoint Species variable"); } } + + if (config->GetKind_Species_Model() == SPECIES_MODEL::FLAMELET) { + + AddVolumeOutput("ADJ_PROGRESS_VARIABLE", "Adjoint_Progress_Variable", "SOLUTION", "Adjoint of the progress variable."); + AddVolumeOutput("ADJ_TOTAL_ENTHALPY", "Adjoint_Total_Enthalpy", "SOLUTION", "Adjoint of the total enthalphy."); + + for (unsigned short i_scalar = 0; i_scalar < config->GetNUserScalars(); i_scalar++) { + const auto& scalar_name = config->GetUserScalarName(i_scalar); + AddVolumeOutput("ADJ_" + scalar_name, "Adjoint_" + scalar_name, "SOLUTION", "Adjoint of " + scalar_name); + } + } } void CAdjFlowOutput::SetVolumeOutputFieldsAdjScalarResidual(const CConfig* config) { @@ -223,12 +296,23 @@ void CAdjFlowOutput::SetVolumeOutputFieldsAdjScalarResidual(const CConfig* confi } } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) { for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { AddVolumeOutput("RES_ADJ_SPECIES_" + std::to_string(iVar), "Residual_Adjoint_Species_" + std::to_string(iVar), "RESIDUAL", "Residual of the adjoint Species variable"); } } + + if (config->GetKind_Species_Model() == SPECIES_MODEL::FLAMELET) { + + AddVolumeOutput("RES_ADJ_PROGRESS_VARIABLE", "Residual_Adjoint_Progress_Variable", "RESIDUAL", "Residual of the adjoint of the progress variable"); + AddVolumeOutput("RES_ADJ_TOTAL_ENTHALPY", "Residual_Adjoint_Total_Enthalpy", "RESIDUAL", "Residual of the adjoint of the total enthalpy"); + + for (unsigned short i_scalar = 0; i_scalar < config->GetNUserScalars(); i_scalar++) { + const auto& scalar_name = config->GetUserScalarName(i_scalar); + AddVolumeOutput("RES_ADJ_" + scalar_name, "Residual_Adjoint_" + scalar_name, "RESIDUAL", "Residual of the adjoint of " + scalar_name); + } + } } void CAdjFlowOutput::LoadVolumeDataAdjScalar(const CConfig* config, const CSolver* const* solver, @@ -258,11 +342,26 @@ void CAdjFlowOutput::LoadVolumeDataAdjScalar(const CConfig* config, const CSolve } } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) { for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { SetVolumeOutputValue("ADJ_SPECIES_" + std::to_string(iVar), iPoint, Node_AdjSpecies->GetSolution(iPoint, iVar)); SetVolumeOutputValue("RES_ADJ_SPECIES_" + std::to_string(iVar), iPoint, Node_AdjSpecies->GetSolution(iPoint, iVar) - Node_AdjSpecies->GetSolution_Old(iPoint, iVar)); } } + + if (config->GetKind_Species_Model() == SPECIES_MODEL::FLAMELET) { + + SetVolumeOutputValue("ADJ_PROGRESS_VARIABLE", iPoint, Node_AdjSpecies->GetSolution(iPoint, I_PROGVAR)); + SetVolumeOutputValue("RES_ADJ_PROGRESS_VARIABLE", iPoint, Node_AdjSpecies->GetSolution(iPoint, I_PROGVAR) - Node_AdjSpecies->GetSolution_Old(iPoint, I_PROGVAR)); + SetVolumeOutputValue("ADJ_TOTAL_ENTHALPY", iPoint, Node_AdjSpecies->GetSolution(iPoint, I_ENTH)); + SetVolumeOutputValue("RES_ADJ_TOTAL_ENTHALPY", iPoint, Node_AdjSpecies->GetSolution(iPoint, I_ENTH) - Node_AdjSpecies->GetSolution_Old(iPoint, I_ENTH)); + + for (unsigned short i_scalar = 0; i_scalar < config->GetNUserScalars(); i_scalar++) { + const auto& scalar_name = config->GetUserScalarName(i_scalar); + SetVolumeOutputValue("ADJ_" + scalar_name, iPoint, Node_AdjSpecies->GetSolution(iPoint, 2 + i_scalar)); + SetVolumeOutputValue("RES_ADJ_" + scalar_name, iPoint, Node_AdjSpecies->GetSolution(iPoint, 2 + i_scalar) - Node_AdjSpecies->GetSolution_Old(iPoint, 2 + i_scalar)); + } + } + } diff --git a/SU2_CFD/src/output/CFlowCompOutput.cpp b/SU2_CFD/src/output/CFlowCompOutput.cpp index a943f4dfb98..4cfe8b9f64c 100644 --- a/SU2_CFD/src/output/CFlowCompOutput.cpp +++ b/SU2_CFD/src/output/CFlowCompOutput.cpp @@ -250,6 +250,8 @@ void CFlowCompOutput::SetVolumeOutputFields(CConfig *config){ AddVolumeOutput("Y_PLUS", "Y_Plus", "PRIMITIVE", "Non-dim. wall distance (Y-Plus)"); } + SetVolumeOutputFieldsScalarPrimitive(config); + //Residuals AddVolumeOutput("RES_DENSITY", "Residual_Density", "RESIDUAL", "Residual of the density"); AddVolumeOutput("RES_MOMENTUM-X", "Residual_Momentum_x", "RESIDUAL", "Residual of the x-momentum component"); @@ -271,6 +273,12 @@ void CFlowCompOutput::SetVolumeOutputFields(CConfig *config){ AddVolumeOutput("LIMITER_ENTHALPY", "Limiter_Enthalpy", "LIMITER", "Limiter value of the enthalpy"); } + SetVolumeOutputFieldsScalarLimiter(config); + + SetVolumeOutputFieldsScalarSource(config); + + SetVolumeOutputFieldsScalarLookup(config); + SetVolumeOutputFieldsScalarMisc(config); // Roe Low Dissipation diff --git a/SU2_CFD/src/output/CFlowIncOutput.cpp b/SU2_CFD/src/output/CFlowIncOutput.cpp index dbf63568bd7..0a500f187ba 100644 --- a/SU2_CFD/src/output/CFlowIncOutput.cpp +++ b/SU2_CFD/src/output/CFlowIncOutput.cpp @@ -38,6 +38,7 @@ CFlowIncOutput::CFlowIncOutput(CConfig *config, unsigned short nDim) : CFlowOutp heat = config->GetEnergy_Equation(); weakly_coupled_heat = config->GetWeakly_Coupled_Heat(); + flamelet = (config->GetKind_Species_Model() == SPECIES_MODEL::FLAMELET); streamwisePeriodic = (config->GetKind_Streamwise_Periodic() != ENUM_STREAMWISE_PERIODIC::NONE); streamwisePeriodic_temperature = config->GetStreamwise_Periodic_Temperature(); @@ -295,7 +296,7 @@ void CFlowIncOutput::SetVolumeOutputFields(CConfig *config){ AddVolumeOutput("VELOCITY-Y", "Velocity_y", "SOLUTION", "y-component of the velocity vector"); if (nDim == 3) AddVolumeOutput("VELOCITY-Z", "Velocity_z", "SOLUTION", "z-component of the velocity vector"); - if (heat || weakly_coupled_heat) + if (heat || weakly_coupled_heat || flamelet) AddVolumeOutput("TEMPERATURE", "Temperature","SOLUTION", "Temperature"); SetVolumeOutputFieldsScalarSolution(config); @@ -331,6 +332,8 @@ void CFlowIncOutput::SetVolumeOutputFields(CConfig *config){ } + SetVolumeOutputFieldsScalarPrimitive(config); + if (config->GetSAParsedOptions().bc) { AddVolumeOutput("INTERMITTENCY", "gamma_BC", "INTERMITTENCY", "Intermittency"); } @@ -352,10 +355,16 @@ void CFlowIncOutput::SetVolumeOutputFields(CConfig *config){ AddVolumeOutput("LIMITER_VELOCITY-Y", "Limiter_Velocity_y", "LIMITER", "Limiter value of the y-velocity"); if (nDim == 3) AddVolumeOutput("LIMITER_VELOCITY-Z", "Limiter_Velocity_z", "LIMITER", "Limiter value of the z-velocity"); - if (heat || weakly_coupled_heat) + if (heat || weakly_coupled_heat || flamelet) AddVolumeOutput("LIMITER_TEMPERATURE", "Limiter_Temperature", "LIMITER", "Limiter value of the temperature"); } + SetVolumeOutputFieldsScalarLimiter(config); + + SetVolumeOutputFieldsScalarSource(config); + + SetVolumeOutputFieldsScalarLookup(config); + SetVolumeOutputFieldsScalarMisc(config); // Streamwise Periodicity @@ -391,7 +400,7 @@ void CFlowIncOutput::LoadVolumeData(CConfig *config, CGeometry *geometry, CSolve if (nDim == 3) SetVolumeOutputValue("VELOCITY-Z", iPoint, Node_Flow->GetSolution(iPoint, 3)); - if (heat) SetVolumeOutputValue("TEMPERATURE", iPoint, Node_Flow->GetSolution(iPoint, nDim+1)); + if (heat || flamelet) SetVolumeOutputValue("TEMPERATURE", iPoint, Node_Flow->GetSolution(iPoint, nDim+1)); if (weakly_coupled_heat) SetVolumeOutputValue("TEMPERATURE", iPoint, Node_Heat->GetSolution(iPoint, 0)); // Radiation solver diff --git a/SU2_CFD/src/output/CFlowOutput.cpp b/SU2_CFD/src/output/CFlowOutput.cpp index 12189f04d69..4f7f94f0d02 100644 --- a/SU2_CFD/src/output/CFlowOutput.cpp +++ b/SU2_CFD/src/output/CFlowOutput.cpp @@ -27,6 +27,8 @@ #include #include +#include +#include #include "../../include/output/CFlowOutput.hpp" @@ -77,7 +79,7 @@ void CFlowOutput::AddAnalyzeSurfaceOutput(const CConfig *config){ } else if (rank == MASTER_NODE) { cout << "\nWARNING: SURFACE_PRESSURE_DROP can only be computed for at least 2 surfaces (outlet, inlet, ...)\n" << endl; } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) { /// DESCRIPTION: Average Species for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { AddHistoryOutput("SURFACE_SPECIES_" + std::to_string(iVar), "Avg_Species_" + std::to_string(iVar), ScreenOutputFormat::FIXED, "SPECIES_COEFF", "Total average species " + std::to_string(iVar) + " on all markers set in MARKER_ANALYZE", HistoryFieldType::COEFFICIENT); @@ -119,7 +121,7 @@ void CFlowOutput::AddAnalyzeSurfaceOutput(const CConfig *config){ AddHistoryOutputPerSurface("SURFACE_TOTAL_TEMPERATURE","Avg_TotalTemp", ScreenOutputFormat::SCIENTIFIC, "FLOW_COEFF_SURF", Marker_Analyze, HistoryFieldType::COEFFICIENT); /// DESCRIPTION: Average total pressure AddHistoryOutputPerSurface("SURFACE_TOTAL_PRESSURE", "Avg_TotalPress", ScreenOutputFormat::SCIENTIFIC, "FLOW_COEFF_SURF", Marker_Analyze, HistoryFieldType::COEFFICIENT); - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) { /// DESCRIPTION: Average Species for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { AddHistoryOutputPerSurface("SURFACE_SPECIES_" + std::to_string(iVar), "Avg_Species_" + std::to_string(iVar), ScreenOutputFormat::FIXED, "SPECIES_COEFF_SURF", Marker_Analyze, HistoryFieldType::COEFFICIENT); @@ -127,7 +129,6 @@ void CFlowOutput::AddAnalyzeSurfaceOutput(const CConfig *config){ /// DESCRIPTION: Species Variance AddHistoryOutputPerSurface("SURFACE_SPECIES_VARIANCE", "Species_Variance", ScreenOutputFormat::SCIENTIFIC, "SPECIES_COEFF_SURF", Marker_Analyze, HistoryFieldType::COEFFICIENT); } - /// END_GROUP } // clang-format on @@ -149,7 +150,7 @@ void CFlowOutput::SetAnalyzeSurface(const CSolver* const*solver, const CGeometry const bool incompressible = config->GetKind_Regime() == ENUM_REGIME::INCOMPRESSIBLE; const bool energy = config->GetEnergy_Equation(); const bool streamwisePeriodic = (config->GetKind_Streamwise_Periodic() != ENUM_STREAMWISE_PERIODIC::NONE); - const bool species = config->GetKind_Species_Model() != SPECIES_MODEL::NONE; + const bool species = config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT; const auto nSpecies = config->GetnSpecies(); const bool axisymmetric = config->GetAxisymmetric(); @@ -378,7 +379,6 @@ void CFlowOutput::SetAnalyzeSurface(const CSolver* const*solver, const CGeometry Allreduce(Surface_MassFlow_Abs_Local, Surface_MassFlow_Abs_Total); Allreduce_su2activematrix(Surface_Species_Local, Surface_Species_Total); - /*--- Compute the value of Surface_Area_Total, and Surface_Pressure_Total, and set the value in the config structure for future use ---*/ @@ -639,7 +639,7 @@ void CFlowOutput::SetAnalyzeSurfaceSpeciesVariance(const CSolver* const*solver, const unsigned short nMarker = config->GetnMarker_All(); const unsigned short Kind_Average = config->GetKind_Average(); - const bool species = config->GetKind_Species_Model() != SPECIES_MODEL::NONE; + const bool species = config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT; const auto nSpecies = config->GetnSpecies(); const bool axisymmetric = config->GetAxisymmetric(); @@ -753,7 +753,7 @@ void CFlowOutput::ConvertVariableSymbolsToIndices(const CPrimitiveIndices& nameToIndex, const std::string& var) { /*--- Primitives of the flow solver. ---*/ @@ -766,10 +766,11 @@ void CFlowOutput::ConvertVariableSymbolsToIndices(const CPrimitiveIndicesGetKind_Species_Model() != SPECIES_MODEL::NONE) { - for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { - AddHistoryOutput("RMS_SPECIES_" + std::to_string(iVar), "rms[rho*Y_" + std::to_string(iVar)+"]", ScreenOutputFormat::FIXED, "RMS_RES", "Root-mean square residual of transported species.", HistoryFieldType::RESIDUAL); + switch (config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: { + for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { + AddHistoryOutput("RMS_SPECIES_" + std::to_string(iVar), "rms[rho*Y_" + std::to_string(iVar)+"]", ScreenOutputFormat::FIXED, "RMS_RES", "Root-mean square residual of transported species.", HistoryFieldType::RESIDUAL); + } + break; } + case SPECIES_MODEL::FLAMELET: { + AddHistoryOutput("RMS_PROGRESS_VARIABLE", "rms[PV]", ScreenOutputFormat::FIXED, "RMS_RES", "Root-mean square residual of the progress variable equation.", HistoryFieldType::RESIDUAL); + AddHistoryOutput("RMS_TOTAL_ENTHALPY", "rms[Enth]", ScreenOutputFormat::FIXED, "RMS_RES", "Root-mean square residual of the total enthalpy equation.", HistoryFieldType::RESIDUAL); + /*--- auxiliary species transport ---*/ + for (auto i_scalar = 0u; i_scalar < config->GetNUserScalars(); i_scalar++){ + const auto& scalar_name = config->GetUserScalarName(i_scalar); + AddHistoryOutput("RMS_"+scalar_name, "rms["+scalar_name+"]", ScreenOutputFormat::FIXED , "RMS_RES", "Root-mean squared residual of the "+scalar_name+" mass fraction equation." , HistoryFieldType::RESIDUAL); + } + break; + } + case SPECIES_MODEL::NONE: break; } } @@ -1022,10 +1037,24 @@ void CFlowOutput::AddHistoryOutputFields_ScalarMAX_RES(const CConfig* config) { break; } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { - for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { - AddHistoryOutput("MAX_SPECIES_" + std::to_string(iVar), "max[rho*Y_" + std::to_string(iVar)+"]", ScreenOutputFormat::FIXED, "MAX_RES", "Maximum residual of transported species.", HistoryFieldType::RESIDUAL); + switch (config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: { + for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { + AddHistoryOutput("MAX_SPECIES_" + std::to_string(iVar), "max[rho*Y_" + std::to_string(iVar)+"]", ScreenOutputFormat::FIXED, "MAX_RES", "Maximum residual of transported species.", HistoryFieldType::RESIDUAL); + } + break; } + case SPECIES_MODEL::FLAMELET: { + AddHistoryOutput("MAX_PROGRESS_VARIABLE", "max[PV]", ScreenOutputFormat::FIXED, "MAX_RES", "Maximum residual of the progress variable equation.", HistoryFieldType::RESIDUAL); + AddHistoryOutput("MAX_TOTAL_ENTHALPY", "max[Enth]", ScreenOutputFormat::FIXED, "MAX_RES", "Maximum residual of the total enthalpy equation.", HistoryFieldType::RESIDUAL); + /*--- auxiliary species transport ---*/ + for (auto i_scalar = 0u; i_scalar < config->GetNUserScalars(); i_scalar++){ + const auto& scalar_name = config->GetUserScalarName(i_scalar); + AddHistoryOutput("MAX_" + scalar_name, "max[" + scalar_name + "]", ScreenOutputFormat::FIXED , "MAX_RES", "Maximum residual of the " + scalar_name + " mass fraction equation." , HistoryFieldType::RESIDUAL); + } + break; + } + case SPECIES_MODEL::NONE: break; } } @@ -1059,10 +1088,24 @@ void CFlowOutput::AddHistoryOutputFields_ScalarBGS_RES(const CConfig* config) { case TURB_TRANS_MODEL::NONE: break; } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { - for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { - AddHistoryOutput("BGS_SPECIES_" + std::to_string(iVar), "bgs[rho*Y_" + std::to_string(iVar)+"]", ScreenOutputFormat::FIXED, "BGS_RES", "BGS residual of transported species.", HistoryFieldType::RESIDUAL); + switch (config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: { + for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { + AddHistoryOutput("BGS_SPECIES_" + std::to_string(iVar), "bgs[rho*Y_" + std::to_string(iVar)+"]", ScreenOutputFormat::FIXED, "BGS_RES", "Maximum residual of transported species.", HistoryFieldType::RESIDUAL); + } + break; } + case SPECIES_MODEL::FLAMELET: { + AddHistoryOutput("BGS_PROGRESS_VARIABLE", "bgs[PV]", ScreenOutputFormat::FIXED, "BGS_RES", "BGS residual of the progress variable equation.", HistoryFieldType::RESIDUAL); + AddHistoryOutput("BGS_TOTAL_ENTHALPY", "bgs[Enth]", ScreenOutputFormat::FIXED, "BGS_RES", "BGS residual of the total enthalpy equation.", HistoryFieldType::RESIDUAL); + /*--- auxiliary species transport ---*/ + for (auto i_scalar = 0u; i_scalar < config->GetNUserScalars(); i_scalar++){ + const auto& scalar_name = config->GetUserScalarName(i_scalar); + AddHistoryOutput("BGS_"+scalar_name, "bgs["+scalar_name+"]", ScreenOutputFormat::FIXED , "BGS_RES", "BGS residual of the "+scalar_name+" mass fraction equation." , HistoryFieldType::RESIDUAL); + } + break; + } + case SPECIES_MODEL::NONE: break; } } @@ -1077,9 +1120,18 @@ void CFlowOutput::AddHistoryOutputFieldsScalarLinsol(const CConfig* config) { AddHistoryOutput("LINSOL_RESIDUAL_TRANS", "LinSolResTrans", ScreenOutputFormat::FIXED, "LINSOL", "Residual of the linear solver for transition solver."); } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { - AddHistoryOutput("LINSOL_ITER_SPECIES", "LinSolIterSpecies", ScreenOutputFormat::INTEGER, "LINSOL", "Number of iterations of the linear solver for species solver."); - AddHistoryOutput("LINSOL_RESIDUAL_SPECIES", "LinSolResSpecies", ScreenOutputFormat::FIXED, "LINSOL", "Residual of the linear solver for species solver."); + switch (config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: { + AddHistoryOutput("LINSOL_ITER_SPECIES", "LinSolIterSpecies", ScreenOutputFormat::INTEGER, "LINSOL", "Number of iterations of the linear solver for species solver."); + AddHistoryOutput("LINSOL_RESIDUAL_SPECIES", "LinSolResSpecies", ScreenOutputFormat::FIXED, "LINSOL", "Residual of the linear solver for species solver."); + break; + } + case SPECIES_MODEL::FLAMELET: { + AddHistoryOutput("LINSOL_ITER_FLAMELET", "LinSolIterSpecies", ScreenOutputFormat::INTEGER, "LINSOL", "Number of iterations of the linear solver for flamelet solver."); + AddHistoryOutput("LINSOL_RESIDUAL_FLAMELET", "LinSolResSpecies", ScreenOutputFormat::FIXED, "LINSOL", "Residual of the linear solver for flamelet solver."); + break; + } + case SPECIES_MODEL::NONE: break; } } // clang-format on @@ -1130,17 +1182,41 @@ void CFlowOutput::LoadHistoryDataScalar(const CConfig* config, const CSolver* co case TURB_TRANS_MODEL::NONE: break; } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { - for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { - SetHistoryOutputValue("RMS_SPECIES_" + std::to_string(iVar), log10(solver[SPECIES_SOL]->GetRes_RMS(iVar))); - SetHistoryOutputValue("MAX_SPECIES_" + std::to_string(iVar), log10(solver[SPECIES_SOL]->GetRes_Max(iVar))); - if (multiZone) { - SetHistoryOutputValue("BGS_SPECIES_" + std::to_string(iVar), log10(solver[SPECIES_SOL]->GetRes_BGS(iVar))); + switch(config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: { + for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { + SetHistoryOutputValue("RMS_SPECIES_" + std::to_string(iVar), log10(solver[SPECIES_SOL]->GetRes_RMS(iVar))); + SetHistoryOutputValue("MAX_SPECIES_" + std::to_string(iVar), log10(solver[SPECIES_SOL]->GetRes_Max(iVar))); + if (multiZone) { + SetHistoryOutputValue("BGS_SPECIES_" + std::to_string(iVar), log10(solver[SPECIES_SOL]->GetRes_BGS(iVar))); + } } + SetHistoryOutputValue("LINSOL_ITER_SPECIES", solver[SPECIES_SOL]->GetIterLinSolver()); + SetHistoryOutputValue("LINSOL_RESIDUAL_SPECIES", log10(solver[SPECIES_SOL]->GetResLinSolver())); + break; } - SetHistoryOutputValue("LINSOL_ITER_SPECIES", solver[SPECIES_SOL]->GetIterLinSolver()); - SetHistoryOutputValue("LINSOL_RESIDUAL_SPECIES", log10(solver[SPECIES_SOL]->GetResLinSolver())); + case SPECIES_MODEL::FLAMELET: { + SetHistoryOutputValue("RMS_PROGRESS_VARIABLE", log10(solver[SPECIES_SOL]->GetRes_RMS(I_PROGVAR))); + SetHistoryOutputValue("MAX_PROGRESS_VARIABLE", log10(solver[SPECIES_SOL]->GetRes_Max(I_PROGVAR))); + SetHistoryOutputValue("RMS_TOTAL_ENTHALPY", log10(solver[SPECIES_SOL]->GetRes_RMS(I_ENTH))); + SetHistoryOutputValue("MAX_TOTAL_ENTHALPY", log10(solver[SPECIES_SOL]->GetRes_Max(I_ENTH))); + /*--- auxiliary species transport ---*/ + for (unsigned short iReactant=0; iReactantGetNUserScalars(); iReactant++){ + const auto& species_name = config->GetUserScalarName(iReactant); + SetHistoryOutputValue("RMS_" + species_name, log10(solver[SPECIES_SOL]->GetRes_RMS(config->GetNControlVars() + iReactant))); + SetHistoryOutputValue("MAX_" + species_name, log10(solver[SPECIES_SOL]->GetRes_Max(config->GetNControlVars() + iReactant))); + if (multiZone) { + SetHistoryOutputValue("BGS_" + species_name, log10(solver[SPECIES_SOL]->GetRes_BGS(config->GetNControlVars() + iReactant))); + } + } + + SetHistoryOutputValue("LINSOL_ITER_FLAMELET", solver[SPECIES_SOL]->GetIterLinSolver()); + SetHistoryOutputValue("LINSOL_RESIDUAL_FLAMELET", log10(solver[SPECIES_SOL]->GetResLinSolver())); + break; + } + + case SPECIES_MODEL::NONE: break; } } @@ -1171,10 +1247,24 @@ void CFlowOutput::SetVolumeOutputFieldsScalarSolution(const CConfig* config){ break; } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { - for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { - AddVolumeOutput("SPECIES_" + std::to_string(iVar), "Species_" + std::to_string(iVar), "SOLUTION", "Species_" + std::to_string(iVar) + " mass fraction"); - } + switch (config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: + for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++){ + AddVolumeOutput("SPECIES_" + std::to_string(iVar), "Species_" + std::to_string(iVar), "SOLUTION", "Species_" + std::to_string(iVar) + " mass fraction"); + } + break; + case SPECIES_MODEL::FLAMELET: + AddVolumeOutput("PROGVAR", "Progress_Variable", "SOLUTION", "Progress variable"); + AddVolumeOutput("ENTHALPY", "Total_Enthalpy", "SOLUTION", "Total enthalpy"); + /*--- auxiliary species ---*/ + for (unsigned short iReactant=0; iReactantGetNUserScalars(); iReactant++) { + const auto& species_name = config->GetUserScalarName(iReactant); + AddVolumeOutput(species_name, species_name, "SOLUTION", species_name + "Mass fraction solution"); + } + + break; + case SPECIES_MODEL::NONE: + break; } } @@ -1195,6 +1285,25 @@ void CFlowOutput::SetVolumeOutputFieldsScalarResidual(const CConfig* config) { break; } + switch (config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: + for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++){ + AddVolumeOutput("RES_SPECIES_" + std::to_string(iVar), "Residual_Species_" + std::to_string(iVar), "RESIDUAL", "Residual of the transported species " + std::to_string(iVar)); + } + break; + case SPECIES_MODEL::FLAMELET: + AddVolumeOutput("RES_PROGVAR", "Residual_Progress_Variable", "RESIDUAL", "Residual of progress variable"); + AddVolumeOutput("RES_ENTHALPY", "Residual_Total_Enthalpy", "RESIDUAL", "Residual of total enthalpy"); + /*--- residuals for auxiliary species transport equations ---*/ + for (unsigned short iReactant=0; iReactantGetNUserScalars(); iReactant++){ + const auto& species_name = config->GetUserScalarName(iReactant); + AddVolumeOutput("RES_" + species_name, "Residual_" + species_name, "RESIDUAL", "Residual of the " + species_name + " equation"); + } + break; + case SPECIES_MODEL::NONE: + break; + } + switch (config->GetKind_Trans_Model()) { case TURB_TRANS_MODEL::LM: AddVolumeOutput("RES_INTERMITTENCY", "Residual_LM_intermittency", "RESIDUAL", "Residual of LM intermittency"); @@ -1204,16 +1313,11 @@ void CFlowOutput::SetVolumeOutputFieldsScalarResidual(const CConfig* config) { case TURB_TRANS_MODEL::NONE: break; } - - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { - for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++){ - AddVolumeOutput("RES_SPECIES_" + std::to_string(iVar), "Residual_Species_" + std::to_string(iVar), "RESIDUAL", "Residual of the transported species " + std::to_string(iVar)); - } - } } -void CFlowOutput::SetVolumeOutputFieldsScalarMisc(const CConfig* config) { - /*--- Place "PRIMITIVE", "LIMITER", and other groups here. ---*/ + +void CFlowOutput::SetVolumeOutputFieldsScalarLimiter(const CConfig* config) { + /*--- Only place outputs of the "SOLUTION" group for species transport here. ---*/ if (config->GetKind_SlopeLimit_Turb() != LIMITER::NONE) { switch (TurbModelFamily(config->GetKind_Turb_Model())) { @@ -1231,15 +1335,41 @@ void CFlowOutput::SetVolumeOutputFieldsScalarMisc(const CConfig* config) { } } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { - if (config->GetKind_SlopeLimit_Species() != LIMITER::NONE) { - for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) - AddVolumeOutput("LIMITER_SPECIES_" + std::to_string(iVar), "Limiter_Species_" + std::to_string(iVar), "LIMITER", "Limiter value of the transported species " + std::to_string(iVar)); - } - for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { - AddVolumeOutput("DIFFUSIVITY_" + std::to_string(iVar), "Diffusivity_" + std::to_string(iVar), "PRIMITIVE", "Diffusivity of the transported species " + std::to_string(iVar)); + if (config->GetKind_SlopeLimit_Species() != LIMITER::NONE) { + switch (config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: + for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) + AddVolumeOutput("LIMITER_SPECIES_" + std::to_string(iVar), "Limiter_Species_" + std::to_string(iVar), "LIMITER", "Limiter value of the transported species " + std::to_string(iVar)); + break; + case SPECIES_MODEL::FLAMELET: + AddVolumeOutput("LIMITER_PROGVAR", "Limiter_Progress_Variable", "LIMITER", "Limiter of progress variable"); + AddVolumeOutput("LIMITER_ENTHALPY", "Limiter_Total_Enthalpy", "LIMITER", "Limiter of total enthalpy"); + /*--- limiter for auxiliary species transport ---*/ + for (unsigned short iReactant=0; iReactant < config->GetNUserScalars(); iReactant++) { + const auto& species_name = config->GetUserScalarName(iReactant); + AddVolumeOutput("LIMITER_" + species_name, "LIMITER_" + species_name, "LIMITER", "Limiter value for the " + species_name + " equation"); + } + break; + case SPECIES_MODEL::NONE: + break; } } +} + +void CFlowOutput::SetVolumeOutputFieldsScalarPrimitive(const CConfig* config) { + /*--- Only place outputs of the "PRIMITIVE" group for scalar transport here. ---*/ + + switch (config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: + for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++){ + AddVolumeOutput("DIFFUSIVITY_" + std::to_string(iVar), "Diffusivity_" + std::to_string(iVar), "PRIMITIVE", "Diffusivity of the transported species " + std::to_string(iVar)); + } + break; + case SPECIES_MODEL::FLAMELET: + break; + case SPECIES_MODEL::NONE: + break; + } switch (config->GetKind_Trans_Model()) { case TURB_TRANS_MODEL::LM: @@ -1256,6 +1386,51 @@ void CFlowOutput::SetVolumeOutputFieldsScalarMisc(const CConfig* config) { AddVolumeOutput("EDDY_VISCOSITY", "Eddy_Viscosity", "PRIMITIVE", "Turbulent eddy viscosity"); } +} + +void CFlowOutput::SetVolumeOutputFieldsScalarSource(const CConfig* config) { + /*--- Only place outputs of the "SOURCE" group for scalar transport here. ---*/ + + switch (config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: + break; + case SPECIES_MODEL::FLAMELET: + AddVolumeOutput("SOURCE_PROGVAR", "Source_Progress_Variable", "SOURCE", "Source Progress Variable"); + /*--- no source term for enthalpy ---*/ + /*--- auxiliary species source terms ---*/ + for (unsigned short iReactant=0; iReactantGetNUserScalars(); iReactant++) { + const auto& species_name = config->GetUserScalarName(iReactant); + AddVolumeOutput("SOURCE_" + species_name, "Source_" + species_name, "SOURCE", "Source " + species_name); + } + break; + case SPECIES_MODEL::NONE: + break; + } +} + + +void CFlowOutput::SetVolumeOutputFieldsScalarLookup(const CConfig* config) { + /*--- Only place outputs of the "LOOKUP" group for scalar transport here. ---*/ + + switch (config->GetKind_Species_Model()) { + case SPECIES_MODEL::SPECIES_TRANSPORT: + break; + case SPECIES_MODEL::FLAMELET: + for (int i_lookup = 0; i_lookup < config->GetNLookups(); ++i_lookup) { + string strname1 = "lookup_" + config->GetLUTLookupName(i_lookup); + AddVolumeOutput(config->GetLUTLookupName(i_lookup), strname1,"LOOKUP", config->GetLUTLookupName(i_lookup)); + } + + AddVolumeOutput("TABLE_MISSES" , "Table_misses" , "LOOKUP", "Lookup table misses"); + break; + case SPECIES_MODEL::NONE: + break; + } +} + +void CFlowOutput::SetVolumeOutputFieldsScalarMisc(const CConfig* config) { + /*--- Only place outputs of the group for scalar transport here that do not fit in other categories. ---*/ + if (config->GetSAParsedOptions().bc) { AddVolumeOutput("INTERMITTENCY", "gamma_BC", "INTERMITTENCY", "Intermittency"); } @@ -1351,16 +1526,57 @@ void CFlowOutput::LoadVolumeDataScalar(const CConfig* config, const CSolver* con SetVolumeOutputValue("WALL_DISTANCE", iPoint, Node_Geo->GetWall_Distance(iPoint)); } - if (config->GetKind_Species_Model() != SPECIES_MODEL::NONE) { - const auto Node_Species = solver[SPECIES_SOL]->GetNodes(); + switch (config->GetKind_Species_Model()) { - for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { - SetVolumeOutputValue("SPECIES_" + std::to_string(iVar), iPoint, Node_Species->GetSolution(iPoint, iVar)); - SetVolumeOutputValue("RES_SPECIES_" + std::to_string(iVar), iPoint, solver[SPECIES_SOL]->LinSysRes(iPoint, iVar)); - SetVolumeOutputValue("DIFFUSIVITY_"+ std::to_string(iVar), iPoint, Node_Species->GetDiffusivity(iPoint,iVar)); - if (config->GetKind_SlopeLimit_Species() != LIMITER::NONE) - SetVolumeOutputValue("LIMITER_SPECIES_" + std::to_string(iVar), iPoint, Node_Species->GetLimiter(iPoint, iVar)); + case SPECIES_MODEL::SPECIES_TRANSPORT: { + const auto Node_Species = solver[SPECIES_SOL]->GetNodes(); + for (unsigned short iVar = 0; iVar < config->GetnSpecies(); iVar++) { + SetVolumeOutputValue("SPECIES_" + std::to_string(iVar), iPoint, Node_Species->GetSolution(iPoint, iVar)); + SetVolumeOutputValue("RES_SPECIES_" + std::to_string(iVar), iPoint, solver[SPECIES_SOL]->LinSysRes(iPoint, iVar)); + SetVolumeOutputValue("DIFFUSIVITY_"+ std::to_string(iVar), iPoint, Node_Species->GetDiffusivity(iPoint,iVar)); + if (config->GetKind_SlopeLimit_Species() != LIMITER::NONE) + SetVolumeOutputValue("LIMITER_SPECIES_" + std::to_string(iVar), iPoint, Node_Species->GetLimiter(iPoint, iVar)); + } + break; + } + + case SPECIES_MODEL::FLAMELET: { + const auto Node_Species = solver[SPECIES_SOL]->GetNodes(); + + SetVolumeOutputValue("PROGVAR", iPoint, Node_Species->GetSolution(iPoint, I_PROGVAR)); + SetVolumeOutputValue("ENTHALPY", iPoint, Node_Species->GetSolution(iPoint, I_ENTH)); + SetVolumeOutputValue("SOURCE_PROGVAR", iPoint, Node_Species->GetScalarSources(iPoint)[I_PROGVAR]); + SetVolumeOutputValue("RES_PROGVAR", iPoint, solver[SPECIES_SOL]->LinSysRes(iPoint, I_PROGVAR)); + SetVolumeOutputValue("RES_ENTHALPY", iPoint, solver[SPECIES_SOL]->LinSysRes(iPoint, I_ENTH)); + SetVolumeOutputValue("TABLE_MISSES", iPoint, su2double(Node_Species->GetTableMisses(iPoint))); + + /*--- auxiliary species transport equations ---*/ + for (unsigned short i_scalar=0; i_scalarGetNUserScalars(); i_scalar++) { + const auto& scalar_name = config->GetUserScalarName(i_scalar); + SetVolumeOutputValue(scalar_name, iPoint, Node_Species->GetSolution(iPoint, config->GetNControlVars() + i_scalar)); + SetVolumeOutputValue("SOURCE_" + scalar_name, iPoint, Node_Species->GetScalarSources(iPoint)[config->GetNControlVars() + i_scalar]); + SetVolumeOutputValue("RES_" + scalar_name, iPoint, solver[SPECIES_SOL]->LinSysRes(iPoint, config->GetNControlVars() + i_scalar)); + } + + if (config->GetKind_SlopeLimit_Species() != LIMITER::NONE) { + SetVolumeOutputValue("LIMITER_PROGVAR", iPoint, Node_Species->GetLimiter(iPoint, I_PROGVAR)); + SetVolumeOutputValue("LIMITER_ENTHALPY", iPoint, Node_Species->GetLimiter(iPoint, I_ENTH)); + /*--- limiter for auxiliary species transport equations ---*/ + for (unsigned short i_scalar=0; i_scalarGetNUserScalars(); i_scalar++) { + const auto& scalar_name = config->GetUserScalarName(i_scalar); + SetVolumeOutputValue("LIMITER_" + scalar_name, iPoint, Node_Species->GetLimiter(iPoint, config->GetNControlVars() + i_scalar)); + } + } + + /*--- variables that we look up from the LUT ---*/ + for (int i_lookup = 0; i_lookup < config->GetNLookups(); ++i_lookup) { + if (config->GetLUTLookupName(i_lookup)!="NULL") + SetVolumeOutputValue(config->GetLUTLookupName(i_lookup), iPoint, Node_Species->GetScalarLookups(iPoint)[i_lookup]); + } + + break; } + case SPECIES_MODEL::NONE: break; } } @@ -2514,6 +2730,10 @@ void CFlowOutput::WriteForcesBreakdown(const CConfig* config, const CSolver* flo << config->GetTemperature_Critical() / config->GetTemperature_Ref() << "\n"; break; + case FLUID_FLAMELET: + file << "Fluid Model: FLAMELET \n"; + break; + case COOLPROP: { CCoolProp auxFluidModel(config->GetFluid_Name()); file << "Fluid Model: CoolProp library \n"; @@ -2833,6 +3053,12 @@ void CFlowOutput::WriteForcesBreakdown(const CConfig* config, const CSolver* flo else file << " psf.\n"; break; + case FLUID_FLAMELET: + file << "Fluid model: FLUID_FLAMELET \n"; + if (si_units) file << " Pa.\n"; + else file << " psf.\n"; + break; + case INC_IDEAL_GAS_POLY: file << "Fluid Model: INC_IDEAL_GAS_POLY \n"; file << "Variable density incompressible flow using ideal gas law.\n"; @@ -2867,6 +3093,13 @@ void CFlowOutput::WriteForcesBreakdown(const CConfig* config, const CSolver* flo file << "Laminar Viscosity (non-dim): " << config->GetMu_ConstantND() << "\n"; break; + case VISCOSITYMODEL::FLAMELET: + file << "Viscosity Model: FLAMELET \n"; + if (si_units) file << " N.s/m^2.\n"; + else file << " lbf.s/ft^2.\n"; + file << "Laminar Viscosity (non-dim): " << config->GetMu_ConstantND() << "\n"; + break; + case VISCOSITYMODEL::COOLPROP: file << "Viscosity Model: CoolProp \n"; break; @@ -2920,6 +3153,12 @@ void CFlowOutput::WriteForcesBreakdown(const CConfig* config, const CSolver* flo file << "Conductivity Model: COOLPROP \n"; break; + case CONDUCTIVITYMODEL::FLAMELET: + file << "Conductivity Model: FLAMELET \n"; + file << "Molecular Conductivity units: " << " W/m^2.K.\n"; + file << "Molecular Conductivity (non-dim): " << config->GetThermal_Conductivity_ConstantND() << "\n"; + break; + case CONDUCTIVITYMODEL::POLYNOMIAL: file << "Viscosity Model: POLYNOMIAL \n"; file << "Kt(T) polynomial coefficients: \n ("; diff --git a/SU2_CFD/src/output/CMultizoneOutput.cpp b/SU2_CFD/src/output/CMultizoneOutput.cpp index a9e0833c422..2c63c04ed2b 100644 --- a/SU2_CFD/src/output/CMultizoneOutput.cpp +++ b/SU2_CFD/src/output/CMultizoneOutput.cpp @@ -173,7 +173,7 @@ void CMultizoneOutput::SetMultizoneHistoryOutputFields(const COutput* const* out /*--- Determine whether Maker_Analyze/Monitoring has to be used. ---*/ auto* Marker = &Marker_Monitoring; - if ((group == "FLOW_COEFF_SURF") || (group == "SPECIES_COEFF_SURF")) + if ((group == "FLOW_COEFF_SURF") || (group == "SPECIES_COEFF_SURF") ) Marker = &Marker_Analyze; else if (group != "AERO_COEFF_SURF" && group != "HEAT_SURF") SU2_MPI::Error("Per Surface output group unknown: " + group, CURRENT_FUNCTION); diff --git a/SU2_CFD/src/output/CNEMOCompOutput.cpp b/SU2_CFD/src/output/CNEMOCompOutput.cpp index d02a2d0c437..f639e3e2d1b 100644 --- a/SU2_CFD/src/output/CNEMOCompOutput.cpp +++ b/SU2_CFD/src/output/CNEMOCompOutput.cpp @@ -261,6 +261,8 @@ void CNEMOCompOutput::SetVolumeOutputFields(CConfig *config){ } + SetVolumeOutputFieldsScalarPrimitive(config); + //Residuals for(iSpecies = 0; iSpecies < nSpecies; iSpecies++) AddVolumeOutput("RES_DENSITY_" + std::to_string(iSpecies), "Residual_Density_" + std::to_string(iSpecies), "RESIDUAL", "Residual of species density " + std::to_string(iSpecies)); @@ -283,6 +285,12 @@ void CNEMOCompOutput::SetVolumeOutputFields(CConfig *config){ AddVolumeOutput("LIMITER_ENERGY", "Limiter_Energy", "LIMITER", "Limiter value of the energy"); } + SetVolumeOutputFieldsScalarLimiter(config); + + SetVolumeOutputFieldsScalarSource(config); + + SetVolumeOutputFieldsScalarLookup(config); + SetVolumeOutputFieldsScalarMisc(config); AddCommonFVMOutputs(config); diff --git a/SU2_CFD/src/output/COutput.cpp b/SU2_CFD/src/output/COutput.cpp index ec61f9754e2..94d7c0a3cb3 100644 --- a/SU2_CFD/src/output/COutput.cpp +++ b/SU2_CFD/src/output/COutput.cpp @@ -263,7 +263,6 @@ void COutput::OutputScreenAndHistory(CConfig *config) { } void COutput::SetupCustomHistoryOutput(const std::string& expression, CustomHistoryOutput& output) const { - std::vector symbols; output.expression = mel::Parse(expression, symbols); diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index 6fe92b813d3..5db213329f9 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -1093,6 +1093,7 @@ void CEulerSolver::SetNondimensionalization(CConfig *config, unsigned short iMes if (viscous) { GetFluidModel()->SetLaminarViscosityModel(config); GetFluidModel()->SetThermalConductivityModel(config); + GetFluidModel()->SetMassDiffusivityModel(config); } } diff --git a/SU2_CFD/src/solvers/CFEM_DG_EulerSolver.cpp b/SU2_CFD/src/solvers/CFEM_DG_EulerSolver.cpp index 877611c3787..7beeb095bbe 100644 --- a/SU2_CFD/src/solvers/CFEM_DG_EulerSolver.cpp +++ b/SU2_CFD/src/solvers/CFEM_DG_EulerSolver.cpp @@ -1092,6 +1092,7 @@ void CFEM_DG_EulerSolver::SetNondimensionalization(CConfig *config, if (viscous) { FluidModel->SetLaminarViscosityModel(config); FluidModel->SetThermalConductivityModel(config); + FluidModel->SetMassDiffusivityModel(config); // nijso: TODO, needs to be tested } if (tkeNeeded) { Energy_FreeStreamND += Tke_FreeStreamND; }; config->SetEnergy_FreeStreamND(Energy_FreeStreamND); diff --git a/SU2_CFD/src/solvers/CIncEulerSolver.cpp b/SU2_CFD/src/solvers/CIncEulerSolver.cpp index e148eb731d3..355b480eb53 100644 --- a/SU2_CFD/src/solvers/CIncEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CIncEulerSolver.cpp @@ -34,6 +34,7 @@ #include "../../include/limiters/CLimiterDetails.hpp" #include "../../../Common/include/toolboxes/geometry_toolbox.hpp" #include "../../include/fluid/CFluidScalar.hpp" +#include "../../include/fluid/CFluidFlamelet.hpp" #include "../../include/fluid/CFluidModel.hpp" @@ -328,6 +329,16 @@ void CIncEulerSolver::SetNondimensionalization(CConfig *config, unsigned short i auxFluidModel->SetTDState_T(Temperature_FreeStream, config->GetSpecies_Init()); break; + case FLUID_FLAMELET: + + config->SetGas_Constant(UNIVERSAL_GAS_CONSTANT / (config->GetMolecular_Weight() / 1000.0)); + Pressure_Thermodynamic = Density_FreeStream * Temperature_FreeStream * config->GetGas_Constant(); + auxFluidModel = new CFluidFlamelet(config, Pressure_Thermodynamic); + auxFluidModel->SetTDState_T(Temperature_FreeStream, config->GetSpecies_Init()); + config->SetPressure_Thermodynamic(Pressure_Thermodynamic); + + break; + default: SU2_MPI::Error("Fluid model not implemented for incompressible solver.", CURRENT_FUNCTION); @@ -403,6 +414,9 @@ void CIncEulerSolver::SetNondimensionalization(CConfig *config, unsigned short i if (tkeNeeded) { Energy_FreeStream += Tke_FreeStream; }; config->SetEnergy_FreeStream(Energy_FreeStream); + /*--- Auxilary (dimensional) FluidModel no longer needed. ---*/ + delete auxFluidModel; + /*--- Compute Mach number ---*/ if (config->GetKind_FluidModel() == CONSTANT_DENSITY) { @@ -449,11 +463,6 @@ void CIncEulerSolver::SetNondimensionalization(CConfig *config, unsigned short i const su2double MassDiffusivityND = config->GetDiffusivity_Constant() / (Velocity_Ref * Length_Ref); config->SetDiffusivity_ConstantND(MassDiffusivityND); - - /*--- Delete the original (dimensional) FluidModel object. No fluid is used for inscompressible cases. ---*/ - - delete auxFluidModel; - /*--- Create one final fluid model object per OpenMP thread to be able to use them in parallel. * GetFluidModel() should be used to automatically access the "right" object of each thread. ---*/ @@ -478,6 +487,11 @@ void CIncEulerSolver::SetNondimensionalization(CConfig *config, unsigned short i fluidModel->SetTDState_T(Temperature_FreeStreamND, config->GetSpecies_Init()); break; + case FLUID_FLAMELET: + fluidModel = new CFluidFlamelet(config, Pressure_Thermodynamic); + fluidModel->SetTDState_T(Temperature_FreeStreamND, config->GetSpecies_Init()); + break; + case INC_IDEAL_GAS_POLY: fluidModel = new CIncIdealGasPolynomial(Gas_ConstantND, Pressure_ThermodynamicND); if (viscous) { @@ -510,7 +524,7 @@ void CIncEulerSolver::SetNondimensionalization(CConfig *config, unsigned short i fluidModel->SetLaminarViscosityModel(config); fluidModel->SetThermalConductivityModel(config); - + fluidModel->SetMassDiffusivityModel(config); } } @@ -625,6 +639,15 @@ void CIncEulerSolver::SetNondimensionalization(CConfig *config, unsigned short i NonDimTable.PrintFooter(); break; + case VISCOSITYMODEL::FLAMELET: + ModelTable << "FLAMELET"; + if (config->GetSystemMeasurements() == SI) Unit << "N.s/m^2"; + else if (config->GetSystemMeasurements() == US) Unit << "lbf.s/ft^2"; + NonDimTable << "Viscosity" << "--" << "--" << Unit.str() << config->GetMu_ConstantND(); + Unit.str(""); + NonDimTable.PrintFooter(); + break; + case VISCOSITYMODEL::COOLPROP: ModelTable << "COOLPROP_VISCOSITY"; if (config->GetSystemMeasurements() == SI) Unit << "N.s/m^2"; @@ -682,6 +705,13 @@ void CIncEulerSolver::SetNondimensionalization(CConfig *config, unsigned short i NonDimTable.PrintFooter(); break; + case CONDUCTIVITYMODEL::FLAMELET: + ModelTable << "FLAMELET"; + Unit << "W/m^2.K"; + NonDimTable << "Molecular Cond." << "--" << "--" << Unit.str() << config->GetThermal_Conductivity_ConstantND(); + Unit.str(""); + break; + case CONDUCTIVITYMODEL::COOLPROP: ModelTable << "COOLPROP"; Unit << "W/m^2.K"; @@ -759,6 +789,23 @@ void CIncEulerSolver::SetNondimensionalization(CConfig *config, unsigned short i NonDimTable.PrintFooter(); break; + case FLUID_FLAMELET: + ModelTable << "FLAMELET"; + Unit << "N.m/kg.K"; + NonDimTable << "Spec. Heat (Cp)" << "--" << "--" << Unit.str() << config->GetSpecific_Heat_CpND(); + Unit.str(""); + Unit << "g/mol"; + NonDimTable << "Molecular weight" << "--" << "--" << Unit.str() << config->GetMolecular_Weight(); + Unit.str(""); + Unit << "N.m/kg.K"; + NonDimTable << "Gas Constant" << "--" << config->GetGas_Constant_Ref() << Unit.str() << config->GetGas_ConstantND(); + Unit.str(""); + Unit << "Pa"; + NonDimTable << "Therm. Pressure" << config->GetPressure_Thermodynamic() << config->GetPressure_Ref() << Unit.str() << config->GetPressure_ThermodynamicND(); + Unit.str(""); + NonDimTable.PrintFooter(); + break; + case INC_IDEAL_GAS_POLY: ModelTable << "INC_IDEAL_GAS_POLY"; Unit.str(""); @@ -2904,8 +2951,6 @@ void CIncEulerSolver::GetOutlet_Properties(CGeometry *geometry, CConfig *config, if (geometry->nodes->GetDomain(iPoint)) { - V_outlet = nodes->GetPrimitive(iPoint); - geometry->vertex[iMarker][iVertex]->GetNormal(Vector); su2double AxiFactor = 1.0; @@ -2921,7 +2966,9 @@ void CIncEulerSolver::GetOutlet_Properties(CGeometry *geometry, CConfig *config, } } - Density = V_outlet[prim_idx.Density()]; + V_outlet = nodes->GetPrimitive(iPoint); + + Density = V_outlet[prim_idx.Density()]; Velocity2 = 0.0; Area = 0.0; MassFlow = 0.0; @@ -2931,6 +2978,7 @@ void CIncEulerSolver::GetOutlet_Properties(CGeometry *geometry, CConfig *config, Velocity2 += Velocity[iDim] * Velocity[iDim]; MassFlow += Vector[iDim] * AxiFactor * Density * Velocity[iDim]; } + Area = sqrt (Area); Outlet_MassFlow[iMarker] += MassFlow; diff --git a/SU2_CFD/src/solvers/CSolver.cpp b/SU2_CFD/src/solvers/CSolver.cpp index 4bd4c0616ee..a7d985811ff 100644 --- a/SU2_CFD/src/solvers/CSolver.cpp +++ b/SU2_CFD/src/solvers/CSolver.cpp @@ -3659,6 +3659,16 @@ void CSolver::LoadInletProfile(CGeometry **geometry, columnValue << config->GetInlet_SpeciesVal(Marker_Tag)[iVar] << "\t"; } break; + case SPECIES_MODEL::FLAMELET: + /*--- 2-equation flamelet model ---*/ + columnName << "PROGRESSVAR" << setw(24) << "ENTHALPYTOT" << setw(24); + columnValue << config->GetInlet_SpeciesVal(Marker_Tag)[0] << "\t" << config->GetInlet_SpeciesVal(Marker_Tag)[1]<<"\t"; + /*--- auxiliary species transport equations ---*/ + for (unsigned short iReactant = 0; iReactant < config->GetNUserScalars(); iReactant++) { + columnName << config->GetUserScalarName(iReactant) << setw(24); + columnValue << config->GetInlet_SpeciesVal(Marker_Tag)[config->GetNControlVars() + iReactant] << "\t"; + } + break; } columnNames.push_back(columnName.str()); diff --git a/SU2_CFD/src/solvers/CSolverFactory.cpp b/SU2_CFD/src/solvers/CSolverFactory.cpp index 02843124a4f..af9669fa392 100644 --- a/SU2_CFD/src/solvers/CSolverFactory.cpp +++ b/SU2_CFD/src/solvers/CSolverFactory.cpp @@ -52,6 +52,7 @@ #include "../../include/solvers/CBaselineSolver_FEM.hpp" #include "../../include/solvers/CRadP1Solver.hpp" #include "../../include/solvers/CSpeciesSolver.hpp" +#include "../../include/solvers/CSpeciesFlameletSolver.hpp" map CSolverFactory::allocatedSolvers; @@ -398,7 +399,10 @@ CSolver* CSolverFactory::CreateSpeciesSolver(CSolver **solver, CGeometry *geomet if (adjoint){ speciesSolver = new CDiscAdjSolver(geometry, config, solver[SPECIES_SOL], RUNTIME_SPECIES_SYS, iMGLevel); } else { - speciesSolver = new CSpeciesSolver(geometry, config, iMGLevel); + if (config->GetKind_Species_Model() == SPECIES_MODEL::SPECIES_TRANSPORT) + speciesSolver = new CSpeciesSolver(geometry, config, iMGLevel); + else if (config->GetKind_Species_Model() == SPECIES_MODEL::FLAMELET) + speciesSolver = new CSpeciesFlameletSolver(geometry, config, iMGLevel); } } return speciesSolver; @@ -411,9 +415,9 @@ CSolver* CSolverFactory::CreateHeatSolver(CSolver **solver, CGeometry *geometry, /*--- Only allocate a heat solver if it should run standalone * or if the weakly coupled heat solver is enabled and no energy equation is included ---*/ - if ((config->GetWeakly_Coupled_Heat() && !config->GetEnergy_Equation()) || config->GetHeatProblem()){ - if (adjoint){ - if (config->GetDiscrete_Adjoint()){ + if ((config->GetWeakly_Coupled_Heat() && !config->GetEnergy_Equation()) || config->GetHeatProblem()) { + if (adjoint) { + if (config->GetDiscrete_Adjoint()) { heatSolver = new CDiscAdjSolver(geometry, config, solver[HEAT_SOL], RUNTIME_HEAT_SYS, iMGLevel); } else { @@ -424,24 +428,24 @@ CSolver* CSolverFactory::CreateHeatSolver(CSolver **solver, CGeometry *geometry, heatSolver = new CHeatSolver(geometry, config, iMGLevel); } } - return heatSolver; + return heatSolver; } CSolver* CSolverFactory::CreateMeshSolver(CSolver **solver, CGeometry *geometry, CConfig *config, int iMGLevel, bool adjoint){ CSolver *meshSolver = nullptr; - if (config->GetDeform_Mesh() && iMGLevel == MESH_0){ - if (!adjoint){ + if (config->GetDeform_Mesh() && iMGLevel == MESH_0) { + if (!adjoint) { meshSolver = new CMeshSolver(geometry, config); } - if (adjoint && config->GetDiscrete_Adjoint()){ + if (adjoint && config->GetDiscrete_Adjoint()) { meshSolver = new CDiscAdjMeshSolver(geometry, config, solver[MESH_SOL]); } } - return meshSolver; + return meshSolver; } CSolver* CSolverFactory::CreateDGSolver(SUB_SOLVER_TYPE kindDGSolver, CGeometry *geometry, CConfig *config, int iMGLevel){ diff --git a/SU2_CFD/src/solvers/CSpeciesFlameletSolver.cpp b/SU2_CFD/src/solvers/CSpeciesFlameletSolver.cpp new file mode 100644 index 00000000000..f79b12b232f --- /dev/null +++ b/SU2_CFD/src/solvers/CSpeciesFlameletSolver.cpp @@ -0,0 +1,511 @@ +/*! + * \file CSpeciesFlameletSolver.cpp + * \brief Main subroutines of CSpeciesFlameletSolver class + * \author D. Mayer, T. Economon, N. Beishuizen + * \version 7.5.1 "Blackbird" + * + * SU2 Project Website: https://su2code.github.io + * + * The SU2 Project is maintained by the SU2 Foundation + * (http://su2foundation.org) + * + * Copyright 2012-2023, SU2 Contributors (cf. AUTHORS.md) + * + * SU2 is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public + * License as published by the Free Software Foundation; either + * version 2.1 of the License, or (at your option) any later version. + * + * SU2 is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with SU2. If not, see . + */ + +#include "../../include/solvers/CSpeciesFlameletSolver.hpp" + +#include "../../../Common/include/parallelization/omp_structure.hpp" +#include "../../../Common/include/toolboxes/geometry_toolbox.hpp" +#include "../../include/fluid/CFluidFlamelet.hpp" +#include "../../include/solvers/CSpeciesSolver.hpp" +#include "../../include/variables/CFlowVariable.hpp" +#include "../../include/variables/CSpeciesFlameletVariable.hpp" + +CSpeciesFlameletSolver::CSpeciesFlameletSolver(CGeometry* geometry, CConfig* config, unsigned short iMesh) + : CSpeciesSolver(geometry, config, true) { + /*--- Dimension of the problem. ---*/ + nVar = config->GetNScalars(); + + Initialize(geometry, config, iMesh, nVar); + + /*--- Initialize the solution to the far-field state everywhere. ---*/ + + nodes = new CSpeciesFlameletVariable(Solution_Inf, nPoint, nDim, nVar, config); + SetBaseClassPointerToNodes(); + + /*--- Store the initial CFL number for all grid points. ---*/ + + const su2double CFL = config->GetCFL(MGLevel) * config->GetCFLRedCoeff_Species(); + SU2_OMP_FOR_STAT(omp_chunk_size) + for (auto iPoint = 0u; iPoint < nPoint; iPoint++) { + nodes->SetLocalCFL(iPoint, CFL); + } + END_SU2_OMP_FOR + Min_CFL_Local = CFL; + Max_CFL_Local = CFL; + Avg_CFL_Local = CFL; + + /*--- Allocates a 3D array with variable "middle" sizes and init to 0. ---*/ + + auto Alloc3D = [](unsigned long M, const vector& N, unsigned long P, vector& X) { + X.resize(M); + for (unsigned long i = 0; i < M; ++i) X[i].resize(N[i], P) = su2double(0.0); + }; + + /*--- Store the values of the temperature and the heat flux density at the boundaries, + used for coupling with a solid donor cell. ---*/ + constexpr auto n_conjugate_var = 4u; + + Alloc3D(nMarker, nVertex, n_conjugate_var, conjugate_var); + for (auto& x : conjugate_var) x = config->GetTemperature_FreeStreamND(); + + /*--- Add the solver name. ---*/ + SolverName = "FLAMELET"; +} + +void CSpeciesFlameletSolver::Preprocessing(CGeometry* geometry, CSolver** solver_container, CConfig* config, + unsigned short iMesh, unsigned short iRKStep, + unsigned short RunTime_EqSystem, bool Output) { + const auto n_user_scalars = config->GetNUserScalars(); + const auto n_control_vars = config->GetNControlVars(); + + vector table_scalar_names(config->GetNScalars()); + table_scalar_names[I_ENTH] = "EnthalpyTot"; + table_scalar_names[I_PROGVAR] = "ProgressVariable"; + + /*--- auxiliary species transport equations---*/ + for (size_t i_aux = 0; i_aux < n_user_scalars; i_aux++) { + table_scalar_names[n_control_vars + i_aux] = config->GetUserScalarName(i_aux); + } + + /*--- we currently only need 1 source term from the LUT for the progress variable + and each auxiliary equations needs 2 source terms ---*/ + unsigned short n_table_sources = 1 + 2 * n_user_scalars; + + vector table_source_names(n_table_sources); + table_source_names[I_SRC_TOT_PROGVAR] = "ProdRateTot_PV"; + /*--- No source term for enthalpy ---*/ + + /*--- For the auxiliary equations, we use a positive (production) and a negative (consumption) term: + S_tot = S_PROD + S_CONS * Y ---*/ + + for (size_t i_aux = 0; i_aux < n_user_scalars; i_aux++) { + /*--- Order of the source terms: S_prod_1, S_cons_1, S_prod_2, S_cons_2, ...---*/ + table_source_names[1 + 2 * i_aux] = config->GetUserSourceName(2 * i_aux); + table_source_names[1 + 2 * i_aux + 1] = config->GetUserSourceName(2 * i_aux + 1); + } + + auto* flowNodes = su2staticcast_p(solver_container[FLOW_SOL]->GetNodes()); + + SU2_OMP_SAFE_GLOBAL_ACCESS(config->SetGlobalParam(config->GetKind_Solver(), RunTime_EqSystem);) + + vector table_lookup_names(config->GetNLookups()); + for (int i_lookup = 0; i_lookup < config->GetNLookups(); ++i_lookup) { + table_lookup_names[i_lookup] = config->GetLUTLookupName(i_lookup); + } + + SU2_OMP_FOR_STAT(omp_chunk_size) + for (auto i_point = 0u; i_point < nPoint; i_point++) { + CFluidModel* fluid_model_local = solver_container[FLOW_SOL]->GetFluidModel(); + su2double* scalars = nodes->GetSolution(i_point); + + /*--- Compute total source terms from the production and consumption. ---*/ + + vector table_sources(n_table_sources); + unsigned long misses = fluid_model_local->GetLookUpTable()->LookUp_XY(table_source_names, table_sources, + scalars[I_PROGVAR], scalars[I_ENTH]); + nodes->SetTableMisses(i_point, misses); + + /*--- The source term for progress variable is always positive, we clip from below to makes sure. --- */ + + vector source_scalar(config->GetNScalars()); + source_scalar[I_PROGVAR] = fmax(EPS, table_sources[I_SRC_TOT_PROGVAR]); + source_scalar[I_ENTH] = 0.0; + + /*--- Source term for the auxiliary species transport equations. ---*/ + for (size_t i_aux = 0; i_aux < config->GetNUserScalars(); i_aux++) { + /*--- The source term for the auxiliary equations consists of a production term and a consumption term: + S_TOT = S_PROD + S_CONS * Y ---*/ + su2double y_aux = scalars[n_control_vars + i_aux]; + su2double source_prod = table_sources[1 + 2 * i_aux]; + su2double source_cons = table_sources[1 + 2 * i_aux + 1]; + source_scalar[n_control_vars + i_aux] = source_prod + source_cons * y_aux; + } + for (auto i_scalar = 0u; i_scalar < nVar; i_scalar++) + nodes->SetScalarSource(i_point, i_scalar, source_scalar[i_scalar]); + + vector lookup_scalar(config->GetNLookups()); + misses = fluid_model_local->GetLookUpTable()->LookUp_XY(table_lookup_names, lookup_scalar, scalars[I_PROGVAR], + scalars[I_ENTH]); + + for (auto i_lookup = 0u; i_lookup < config->GetNLookups(); i_lookup++) { + nodes->SetLookupScalar(i_point, lookup_scalar[i_lookup], i_lookup); + } + + su2double T = flowNodes->GetTemperature(i_point); + fluid_model_local->SetTDState_T(T, scalars); + /*--- set the diffusivity in the fluid model to the diffusivity obtained from the lookup table ---*/ + for (auto i_scalar = 0u; i_scalar < nVar; ++i_scalar) { + nodes->SetDiffusivity(i_point, fluid_model_local->GetMassDiffusivity(i_scalar), i_scalar); + } + + if (!Output) LinSysRes.SetBlock_Zero(i_point); + } + END_SU2_OMP_FOR + + /*--- Clear Residual and Jacobian. Upwind second order reconstruction and gradients ---*/ + CommonPreprocessing(geometry, config, Output); +} + +void CSpeciesFlameletSolver::SetInitialCondition(CGeometry** geometry, CSolver*** solver_container, CConfig* config, + unsigned long ExtIter) { + const bool Restart = (config->GetRestart() || config->GetRestart_Flow()); + + if ((!Restart) && ExtIter == 0) { + if (rank == MASTER_NODE) { + cout << "Initializing progress variable and total enthalpy (using temperature)" << endl; + } + + vector scalar_init(nVar, 0.0); + const su2double* flame_init = config->GetFlameInit(); + const su2double flame_offset[3] = {flame_init[0], flame_init[1], flame_init[2]}; + const su2double flame_normal[3] = {flame_init[3], flame_init[4], flame_init[5]}; + const su2double flame_thickness = flame_init[6]; + const su2double flame_burnt_thickness = flame_init[7]; + + const su2double flamenorm = GeometryToolbox::Norm(nDim, flame_normal); + const su2double temp_inlet = config->GetInc_Temperature_Init(); + su2double prog_inlet = config->GetSpecies_Init()[I_PROGVAR]; + su2double enth_inlet = config->GetSpecies_Init()[I_ENTH]; + + su2double prog_burnt; + su2double prog_unburnt = 0.0; + + if (rank == MASTER_NODE) { + cout << "initial condition: T = " << temp_inlet << endl; + cout << "initial condition: c = " << prog_inlet << endl; + cout << "initial condition: h = " << enth_inlet << endl; + } + + su2double point_loc; + + CFluidModel* fluid_model_local; + + vector look_up_tags; + vector look_up_data; + + unsigned long n_not_iterated_local = 0; + unsigned long n_not_in_domain_local = 0; + unsigned long n_points_unburnt_local = 0; + unsigned long n_points_burnt_local = 0; + unsigned long n_points_flame_local = 0; + unsigned long n_not_iterated_global; + unsigned long n_not_in_domain_global; + unsigned long n_points_burnt_global; + unsigned long n_points_flame_global; + unsigned long n_points_unburnt_global; + + for (unsigned long i_mesh = 0; i_mesh <= config->GetnMGLevels(); i_mesh++) { + fluid_model_local = solver_container[i_mesh][FLOW_SOL]->GetFluidModel(); + + prog_burnt = *fluid_model_local->GetLookUpTable()->GetTableLimitsX().second; + for (unsigned long i_point = 0; i_point < nPointDomain; i_point++) { + auto coords = geometry[i_mesh]->nodes->GetCoord(i_point); + + /*--- Determine if point is above or below the plane, assuming the normal + is pointing towards the burned region. ---*/ + point_loc = 0.0; + for (unsigned short i_dim = 0; i_dim < geometry[i_mesh]->GetnDim(); i_dim++) { + point_loc += flame_normal[i_dim] * (coords[i_dim] - flame_offset[i_dim]); + } + + /*--- Compute the exact distance from point to plane. ---*/ + point_loc = point_loc / flamenorm; + + /* --- Unburnt region upstream of the flame. --- */ + if (point_loc <= 0) { + scalar_init[I_PROGVAR] = prog_unburnt; + n_points_unburnt_local++; + + /* --- Flame zone where we lineary increase from unburnt to burnt conditions. --- */ + } else if ((point_loc > 0) && (point_loc <= flame_thickness)) { + scalar_init[I_PROGVAR] = prog_unburnt + point_loc * (prog_burnt - prog_unburnt) / flame_thickness; + n_points_flame_local++; + + /* --- Burnt region behind the flame zone. --- */ + } else if ((point_loc > flame_thickness) && (point_loc <= flame_thickness + flame_burnt_thickness)) { + scalar_init[I_PROGVAR] = prog_burnt; + n_points_burnt_local++; + + /* --- Unburnt region downstream of the flame and burnt region. --- */ + } else { + scalar_init[I_PROGVAR] = prog_unburnt; + n_points_unburnt_local++; + } + + n_not_iterated_local += fluid_model_local->GetEnthFromTemp(enth_inlet, prog_inlet, temp_inlet, enth_inlet); + scalar_init[I_ENTH] = enth_inlet; + + n_not_in_domain_local += fluid_model_local->GetLookUpTable()->LookUp_XY( + look_up_tags, look_up_data, scalar_init[I_PROGVAR], scalar_init[I_ENTH]); + + /* --- Initialize the auxiliary transported scalars (not controlling variables). --- */ + for (int i_scalar = config->GetNControlVars(); i_scalar < config->GetNScalars(); ++i_scalar) { + scalar_init[i_scalar] = config->GetSpecies_Init()[i_scalar]; + } + + solver_container[i_mesh][SPECIES_SOL]->GetNodes()->SetSolution(i_point, scalar_init.data()); + } + + solver_container[i_mesh][SPECIES_SOL]->InitiateComms(geometry[i_mesh], config, SOLUTION); + solver_container[i_mesh][SPECIES_SOL]->CompleteComms(geometry[i_mesh], config, SOLUTION); + + solver_container[i_mesh][FLOW_SOL]->InitiateComms(geometry[i_mesh], config, SOLUTION); + solver_container[i_mesh][FLOW_SOL]->CompleteComms(geometry[i_mesh], config, SOLUTION); + + solver_container[i_mesh][FLOW_SOL]->Preprocessing(geometry[i_mesh], solver_container[i_mesh], config, i_mesh, + NO_RK_ITER, RUNTIME_FLOW_SYS, false); + } + + /* --- Sum up some global counters over processes. --- */ + SU2_MPI::Reduce(&n_not_in_domain_local, &n_not_in_domain_global, 1, MPI_UNSIGNED_LONG, MPI_SUM, MASTER_NODE, + SU2_MPI::GetComm()); + SU2_MPI::Reduce(&n_not_iterated_local, &n_not_iterated_global, 1, MPI_UNSIGNED_LONG, MPI_SUM, MASTER_NODE, + SU2_MPI::GetComm()); + SU2_MPI::Reduce(&n_points_unburnt_local, &n_points_unburnt_global, 1, MPI_UNSIGNED_LONG, MPI_SUM, MASTER_NODE, + SU2_MPI::GetComm()); + SU2_MPI::Reduce(&n_points_burnt_local, &n_points_burnt_global, 1, MPI_UNSIGNED_LONG, MPI_SUM, MASTER_NODE, + SU2_MPI::GetComm()); + SU2_MPI::Reduce(&n_points_flame_local, &n_points_flame_global, 1, MPI_UNSIGNED_LONG, MPI_SUM, MASTER_NODE, + SU2_MPI::GetComm()); + + if (rank == MASTER_NODE) { + cout << endl; + cout << " Number of points in unburnt region: " << n_points_unburnt_global << "." << endl; + cout << " Number of points in burnt region : " << n_points_burnt_global << "." << endl; + cout << " Number of points in flame zone : " << n_points_flame_global << "." << endl; + + if (n_not_in_domain_global > 0) + cout << " Initial condition: Number of points outside of table domain: " << n_not_in_domain_global << " !!!" + << endl; + + if (n_not_iterated_global > 0) + cout << " Initial condition: Number of points in which enthalpy could not be iterated: " << n_not_iterated_global + << " !!!" << endl; + } + } +} + +void CSpeciesFlameletSolver::SetPreconditioner(CGeometry* geometry, CSolver** solver_container, CConfig* config) { + const bool variable_density = (config->GetKind_DensityModel() == INC_DENSITYMODEL::VARIABLE); + const bool implicit = (config->GetKind_TimeIntScheme_Flow() == EULER_IMPLICIT); + + SU2_OMP_FOR_STAT(omp_chunk_size) + for (auto iPoint = 0u; iPoint < nPointDomain; iPoint++) { + /*--- Access the primitive variables at this node. ---*/ + + su2double Density = solver_container[FLOW_SOL]->GetNodes()->GetDensity(iPoint); + su2double BetaInc2 = solver_container[FLOW_SOL]->GetNodes()->GetBetaInc2(iPoint); + su2double Temperature = solver_container[FLOW_SOL]->GetNodes()->GetTemperature(iPoint); + + su2double SolP = solver_container[FLOW_SOL]->LinSysSol(iPoint, prim_idx.Pressure()); + su2double SolT = solver_container[FLOW_SOL]->LinSysSol(iPoint, prim_idx.Temperature()); + + /*--- We need the derivative of the equation of state to build the + preconditioning matrix. For now, the only option is the ideal gas + law, but in the future, dRhodT should be in the fluid model. ---*/ + + su2double dRhodT = 0.0; + if (variable_density) { + dRhodT = -Density / Temperature; + } + + /*--- Passive scalars have no impact on the density. ---*/ + + su2double dRhodC = 0.0; + + /*--- Modify matrix diagonal with term including volume and time step. ---*/ + + su2double Vol = geometry->nodes->GetVolume(iPoint); + su2double Delta = + Vol / (config->GetCFLRedCoeff_Species() * solver_container[FLOW_SOL]->GetNodes()->GetDelta_Time(iPoint)); + + /*--- Calculating the inverse of the preconditioning matrix + that multiplies the time derivative during time integration. ---*/ + + if (implicit) { + for (unsigned short iVar = 0; iVar < nVar; iVar++) { + su2double scalar = nodes->GetSolution(iPoint, iVar); + + /*--- Compute the lag terms for the decoupled linear system from + the mean flow equations and add to the residual for the scalar. + In short, we are effectively making these terms explicit. ---*/ + + su2double artcompc1 = SolP * scalar / (Density * BetaInc2); + su2double artcompc2 = SolT * dRhodT * scalar / (Density); + + LinSysRes(iPoint, iVar) += artcompc1 + artcompc2; + + /*--- Add the extra Jacobian term to the scalar system. ---*/ + + su2double Jaccomp = scalar * dRhodC + Density; + su2double JacTerm = Jaccomp * Delta; + + Jacobian.AddVal2Diag(iPoint, iVar, JacTerm); + } + } + } + END_SU2_OMP_FOR +} + +void CSpeciesFlameletSolver::Source_Residual(CGeometry* geometry, CSolver** solver_container, + CNumerics** numerics_container, CConfig* config, unsigned short iMesh) { + SU2_OMP_FOR_STAT(omp_chunk_size) + for (auto i_point = 0u; i_point < nPointDomain; i_point++) { + + /*--- Add source terms from the lookup table directly to the residual. ---*/ + + for (auto i_var = 0; i_var < nVar; i_var++) { + LinSysRes(i_point, i_var) -= nodes->GetScalarSources(i_point)[i_var] * geometry->nodes->GetVolume(i_point); + } + } + END_SU2_OMP_FOR + + /*--- call the species solver for the shared sources (axisymmetric) ---*/ + CSpeciesSolver::Source_Residual(geometry, solver_container, numerics_container, config, iMesh); +} + +void CSpeciesFlameletSolver::BC_Inlet(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, + CNumerics* visc_numerics, CConfig* config, unsigned short val_marker) { + string Marker_Tag = config->GetMarker_All_TagBound(val_marker); + + su2double temp_inlet = config->GetInlet_Ttotal(Marker_Tag); + const su2double* scalar_inlet = config->GetInlet_SpeciesVal(Marker_Tag); + + /*--- We compute inlet enthalpy from the temperature and progress variable. ---*/ + su2double enth_inlet = scalar_inlet[I_ENTH]; + solver_container[FLOW_SOL]->GetFluidModel()->GetEnthFromTemp(enth_inlet, scalar_inlet[I_PROGVAR], temp_inlet, + scalar_inlet[I_ENTH]); + + SU2_OMP_FOR_STAT(OMP_MIN_SIZE) + for (auto iVertex = 0u; iVertex < geometry->nVertex[val_marker]; iVertex++) { + Inlet_SpeciesVars[val_marker][iVertex][I_ENTH] = enth_inlet; + END_SU2_OMP_FOR + + } + + /*--- Call the general inlet boundary condition implementation. ---*/ + CSpeciesSolver::BC_Inlet(geometry, solver_container, conv_numerics, visc_numerics, config, val_marker); +} + +void CSpeciesFlameletSolver::BC_Isothermal_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 string Marker_Tag = config->GetMarker_All_TagBound(val_marker); + su2double temp_wall = config->GetIsothermal_Temperature(Marker_Tag); + CFluidModel* fluid_model_local = solver_container[FLOW_SOL]->GetFluidModel(); + auto* flowNodes = su2staticcast_p(solver_container[FLOW_SOL]->GetNodes()); + su2double enth_init, enth_wall, prog_wall; + unsigned long n_not_iterated = 0; + + /*--- Loop over all the vertices on this boundary marker. ---*/ + SU2_OMP_FOR_STAT(OMP_MIN_SIZE) + for (unsigned long iVertex = 0; iVertex < geometry->nVertex[val_marker]; iVertex++) { + unsigned long iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); + + if (cht_mode) temp_wall = GetConjugateHeatVariable(val_marker, iVertex, 0); + + /*--- Check if the node belongs to the domain (i.e., not a halo node). ---*/ + + if (geometry->nodes->GetDomain(iPoint)) { + if (config->GetMarker_StrongBC(Marker_Tag) == true) { + /*--- Initial guess for enthalpy value. ---*/ + + enth_init = nodes->GetSolution(iPoint, I_ENTH); + enth_wall = enth_init; + + /*--- Set enthalpy on the wall. ---*/ + + prog_wall = solver_container[SPECIES_SOL]->GetNodes()->GetSolution(iPoint)[I_PROGVAR]; + n_not_iterated += fluid_model_local->GetEnthFromTemp(enth_wall, prog_wall, temp_wall, enth_init); + + /*--- Impose the value of the enthalpy as a strong boundary + condition (Dirichlet) and remove any + contribution to the residual at this node. ---*/ + + nodes->SetSolution(iPoint, I_ENTH, enth_wall); + nodes->SetSolution_Old(iPoint, I_ENTH, enth_wall); + + LinSysRes(iPoint, I_ENTH) = 0.0; + + nodes->SetVal_ResTruncError_Zero(iPoint, I_ENTH); + + if (implicit) { + unsigned long total_index = iPoint * nVar + I_ENTH; + Jacobian.DeleteValsRowi(total_index); + } + } else { + /*--- Weak BC formulation. ---*/ + const auto Normal = geometry->vertex[val_marker][iVertex]->GetNormal(); + + const su2double Area = GeometryToolbox::Norm(nDim, Normal); + + const auto Point_Normal = geometry->vertex[val_marker][iVertex]->GetNormal_Neighbor(); + + /*--- Get coordinates of i & nearest normal and compute distance. ---*/ + + const auto Coord_i = geometry->nodes->GetCoord(iPoint); + const auto Coord_j = geometry->nodes->GetCoord(Point_Normal); + su2double Edge_Vector[MAXNDIM]; + GeometryToolbox::Distance(nDim, Coord_j, Coord_i, Edge_Vector); + su2double dist_ij_2 = GeometryToolbox::SquaredNorm(nDim, Edge_Vector); + su2double dist_ij = sqrt(dist_ij_2); + + /*--- Compute the normal gradient in temperature using Twall. ---*/ + + su2double dTdn = -(flowNodes->GetTemperature(Point_Normal) - temp_wall) / dist_ij; + + /*--- Get thermal conductivity. ---*/ + + su2double thermal_conductivity = flowNodes->GetThermalConductivity(iPoint); + + /*--- Apply a weak boundary condition for the energy equation. + Compute the residual due to the prescribed heat flux. ---*/ + + LinSysRes(iPoint, I_ENTH) -= thermal_conductivity * dTdn * Area; + } + } + } + END_SU2_OMP_FOR + + if (rank == MASTER_NODE && n_not_iterated > 0) { + cout << " !!! Wall bc (" << Marker_Tag + << "): Number of points in which enthalpy could not be iterated: " << n_not_iterated << " !!!" << endl; + } +} + +void CSpeciesFlameletSolver::BC_Isothermal_Wall(CGeometry* geometry, CSolver** solver_container, + CNumerics* conv_numerics, CNumerics* visc_numerics, CConfig* config, + unsigned short val_marker) { + BC_Isothermal_Wall_Generic(geometry, solver_container, conv_numerics, visc_numerics, config, val_marker); +} + +void CSpeciesFlameletSolver::BC_ConjugateHeat_Interface(CGeometry* geometry, CSolver** solver_container, + CNumerics* conv_numerics, CConfig* config, + unsigned short val_marker) { + BC_Isothermal_Wall_Generic(geometry, solver_container, conv_numerics, nullptr, config, val_marker, true); +} diff --git a/SU2_CFD/src/solvers/CSpeciesSolver.cpp b/SU2_CFD/src/solvers/CSpeciesSolver.cpp index 43a88edf3ea..5acebd6164a 100644 --- a/SU2_CFD/src/solvers/CSpeciesSolver.cpp +++ b/SU2_CFD/src/solvers/CSpeciesSolver.cpp @@ -36,12 +36,64 @@ template class CScalarSolver; CSpeciesSolver::CSpeciesSolver(CGeometry* geometry, CConfig* config, unsigned short iMesh) : CScalarSolver(geometry, config, true) { - /*--- Store if an implicit scheme is used, for use during periodic boundary conditions. ---*/ - SetImplicitPeriodic(config->GetKind_TimeIntScheme_Species() == EULER_IMPLICIT); /*--- Dimension of the problem. ---*/ nVar = config->GetnSpecies(); + + Initialize(geometry, config, iMesh, nVar); + + /*--- Initialize the solution to the far-field state everywhere. ---*/ + + nodes = new CSpeciesVariable(Solution_Inf, nPoint, nDim, nVar, config); + SetBaseClassPointerToNodes(); + + /*--- Initialize the mass diffusivity. Nondimensionalization done in the flow solver. ---*/ + SU2_OMP_FOR_STAT(omp_chunk_size) + for (unsigned long iPoint = 0; iPoint < nPoint; iPoint++) { + for (auto iVar = 0u; iVar <= nVar; iVar++) { + const auto MassDiffusivity = config->GetDiffusivity_ConstantND(); + nodes->SetDiffusivity(iPoint, MassDiffusivity, iVar); + } + } + END_SU2_OMP_FOR + + /*--- MPI solution ---*/ + + InitiateComms(geometry, config, SOLUTION); + CompleteComms(geometry, config, SOLUTION); + + SlidingState.resize(nMarker); + SlidingStateNodes.resize(nMarker); + + for (unsigned long iMarker = 0; iMarker < nMarker; iMarker++) { + if (config->GetMarker_All_KindBC(iMarker) == FLUID_INTERFACE) { + SlidingState[iMarker].resize(nVertex[iMarker], nPrimVar+1) = nullptr; + SlidingStateNodes[iMarker].resize(nVertex[iMarker],0); + } + } + + /*--- Store the initial CFL number for all grid points. ---*/ + + const su2double CFL = config->GetCFL(MGLevel) * config->GetCFLRedCoeff_Species(); + SU2_OMP_FOR_STAT(omp_chunk_size) + for (auto iPoint = 0u; iPoint < nPoint; iPoint++) { + nodes->SetLocalCFL(iPoint, CFL); + } + END_SU2_OMP_FOR + Min_CFL_Local = CFL; + Max_CFL_Local = CFL; + Avg_CFL_Local = CFL; + + /*--- Add the solver name. ---*/ + SolverName = "SPECIES"; +} + + +void CSpeciesSolver::Initialize(CGeometry* geometry, CConfig* config, unsigned short iMesh, unsigned short nVar) { + /*--- Store if an implicit scheme is used, for use during periodic boundary conditions. ---*/ + SetImplicitPeriodic(config->GetKind_TimeIntScheme_Species() == EULER_IMPLICIT); + nPrimVar = nVar; if (nVar > MAXNVAR) @@ -58,9 +110,8 @@ CSpeciesSolver::CSpeciesSolver(CGeometry* geometry, CConfig* config, unsigned sh nDim = geometry->GetnDim(); - /*--- Single grid simulation ---*/ - if (iMesh == MESH_0 || config->GetMGCycle() == FULLMG_CYCLE) { +if (iMesh == MESH_0 || config->GetMGCycle() == FULLMG_CYCLE) { /*--- Define some auxiliary vector related with the residual ---*/ @@ -108,41 +159,11 @@ CSpeciesSolver::CSpeciesSolver(CGeometry* geometry, CConfig* config, unsigned sh Solution_Inf[iVar] = config->GetSpecies_Init()[iVar]; } - /*--- Initialize the solution to the far-field state everywhere. ---*/ - - nodes = new CSpeciesVariable(Solution_Inf, nPoint, nDim, nVar, config); - SetBaseClassPointerToNodes(); - - /*--- Initialize the mass diffusivity. Nondimensionalization done in the flow solver. ---*/ - SU2_OMP_FOR_STAT(omp_chunk_size) - for (unsigned long iPoint = 0; iPoint < nPoint; iPoint++) { - for (auto iVar = 0u; iVar <= nVar; iVar++) { - const auto MassDiffusivity = config->GetDiffusivity_ConstantND(); - nodes->SetDiffusivity(iPoint, MassDiffusivity, iVar); - } - } - END_SU2_OMP_FOR - - /*--- MPI solution ---*/ - - InitiateComms(geometry, config, SOLUTION); - CompleteComms(geometry, config, SOLUTION); - - SlidingState.resize(nMarker); - SlidingStateNodes.resize(nMarker); - - for (unsigned long iMarker = 0; iMarker < nMarker; iMarker++) { - if (config->GetMarker_All_KindBC(iMarker) == FLUID_INTERFACE) { - SlidingState[iMarker].resize(nVertex[iMarker], nPrimVar+1) = nullptr; - SlidingStateNodes[iMarker].resize(nVertex[iMarker],0); - } - } - /*--- Set the column number for species in inlet-files. * e.g. Coords(nDim), Temp(1), VelMag(1), Normal(nDim), Turb(1 or 2), Species(arbitrary) ---*/ Inlet_Position = nDim + 2 + nDim + config->GetnTurbVar(); - /*-- Allocation of inlet-values. Will be filled either by an inlet files, + /*-- Allocation of inlet-values. Will be filled either by an inlet file, * or uniformly by a uniform boundary condition. ---*/ Inlet_SpeciesVars.resize(nMarker); @@ -154,23 +175,9 @@ CSpeciesSolver::CSpeciesSolver(CGeometry* geometry, CConfig* config, unsigned sh } } } - - /*--- Store the initial CFL number for all grid points. ---*/ - - const su2double CFL = config->GetCFL(MGLevel) * config->GetCFLRedCoeff_Species(); - SU2_OMP_FOR_STAT(omp_chunk_size) - for (auto iPoint = 0u; iPoint < nPoint; iPoint++) { - nodes->SetLocalCFL(iPoint, CFL); - } - END_SU2_OMP_FOR - Min_CFL_Local = CFL; - Max_CFL_Local = CFL; - Avg_CFL_Local = CFL; - - /*--- Add the solver name. ---*/ - SolverName = "SPECIES"; } + void CSpeciesSolver::LoadRestart(CGeometry** geometry, CSolver*** solver, CConfig* config, int val_iter, bool val_update_geo) { /*--- Restart the solution from file information ---*/ @@ -199,9 +206,12 @@ void CSpeciesSolver::LoadRestart(CGeometry** geometry, CSolver*** solver, CConfi const bool incompressible = (config->GetKind_Regime() == ENUM_REGIME::INCOMPRESSIBLE); const bool energy = config->GetEnergy_Equation(); + const bool flamelet = (config->GetKind_FluidModel() == FLUID_FLAMELET); const bool weakly_coupled_heat = config->GetWeakly_Coupled_Heat(); - if (incompressible && ((!energy) && (!weakly_coupled_heat))) skipVars--; + /*--- for the flamelet model, the temperature is saved to file, but the energy equation is off ---*/ + + if (incompressible && ((!energy) && (!weakly_coupled_heat) && (!flamelet))) skipVars--; /*--- Load data from the restart into correct containers. ---*/ @@ -322,7 +332,6 @@ void CSpeciesSolver::Viscous_Residual(unsigned long iEdge, CGeometry* geometry, void CSpeciesSolver::BC_Inlet(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, CNumerics* visc_numerics, CConfig* config, unsigned short val_marker) { - string Marker_Tag = config->GetMarker_All_TagBound(val_marker); const bool implicit = (config->GetKind_TimeIntScheme() == EULER_IMPLICIT); @@ -336,7 +345,10 @@ void CSpeciesSolver::BC_Inlet(CGeometry* geometry, CSolver** solver_container, C if (!geometry->nodes->GetDomain(iPoint)) continue; - if (config->GetSpecies_StrongBC()) { + /*--- Identify the boundary by string name ---*/ + string Marker_Tag = config->GetMarker_All_TagBound(val_marker); + + if (config->GetMarker_StrongBC(Marker_Tag)==true) { nodes->SetSolution_Old(iPoint, Inlet_SpeciesVars[val_marker][iVertex]); LinSysRes.SetBlock_Zero(iPoint); @@ -426,7 +438,7 @@ su2double CSpeciesSolver::GetInletAtVertex(su2double *val_inlet, if ((config->GetMarker_All_KindBC(iMarker) == INLET_FLOW) && (config->GetMarker_All_TagBound(iMarker) == val_marker)) { - for (iVertex = 0; iVertex < nVertex[iMarker]; iVertex++){ + for (iVertex = 0; iVertex < nVertex[iMarker]; iVertex++) { iPoint = geometry->vertex[iMarker][iVertex]->GetNode(); @@ -481,14 +493,17 @@ void CSpeciesSolver::BC_Outlet(CGeometry* geometry, CSolver** solver_container, SU2_OMP_FOR_STAT(OMP_MIN_SIZE) for (auto iVertex = 0u; iVertex < geometry->nVertex[val_marker]; iVertex++) { - /* strong zero flux Neumann boundary condition at the outlet */ + /*--- Strong zero flux Neumann boundary condition at the outlet ---*/ const auto iPoint = geometry->vertex[val_marker][iVertex]->GetNode(); /*--- Check if the node belongs to the domain (i.e., not a halo node) ---*/ if (!geometry->nodes->GetDomain(iPoint)) continue; - if (config->GetSpecies_StrongBC()) { + /*--- Identify the boundary by string name ---*/ + string Marker_Tag = config->GetMarker_All_TagBound(val_marker); + + if (config->GetMarker_StrongBC(Marker_Tag)==true) { /*--- Allocate the value at the outlet ---*/ auto Point_Normal = geometry->vertex[val_marker][iVertex]->GetNormal_Neighbor(); diff --git a/SU2_CFD/src/solvers/CTurbSolver.cpp b/SU2_CFD/src/solvers/CTurbSolver.cpp index d68c550d185..ffd7152346b 100644 --- a/SU2_CFD/src/solvers/CTurbSolver.cpp +++ b/SU2_CFD/src/solvers/CTurbSolver.cpp @@ -129,8 +129,9 @@ void CTurbSolver::LoadRestart(CGeometry** geometry, CSolver*** solver, CConfig* const bool incompressible = (config->GetKind_Regime() == ENUM_REGIME::INCOMPRESSIBLE); const bool energy = config->GetEnergy_Equation(); const bool weakly_coupled_heat = config->GetWeakly_Coupled_Heat(); + const bool flamelet = (config->GetKind_FluidModel() == FLUID_FLAMELET); - if (incompressible && ((!energy) && (!weakly_coupled_heat))) skipVars--; + if (incompressible && ((!energy) && (!weakly_coupled_heat) && (!flamelet))) skipVars--; /*--- Load data from the restart into correct containers. ---*/ diff --git a/SU2_CFD/src/variables/CIncNSVariable.cpp b/SU2_CFD/src/variables/CIncNSVariable.cpp index c000575c218..c132cb20efd 100644 --- a/SU2_CFD/src/variables/CIncNSVariable.cpp +++ b/SU2_CFD/src/variables/CIncNSVariable.cpp @@ -58,8 +58,8 @@ bool CIncNSVariable::SetPrimVar(unsigned long iPoint, su2double eddy_visc, su2do /*--- Set the value of the temperature directly ---*/ - su2double Temperature = Solution(iPoint, nDim+1); - const auto check_temp = SetTemperature(iPoint,Temperature); + su2double Temperature = Solution(iPoint, indices.Temperature()); + auto check_temp = SetTemperature(iPoint, Temperature); /*--- Use the fluid model to compute the new value of density. Note that the thermodynamic pressure is constant and decoupled @@ -69,6 +69,12 @@ bool CIncNSVariable::SetPrimVar(unsigned long iPoint, su2double eddy_visc, su2do FluidModel->SetTDState_T(Temperature, scalar); + /*--- for FLAMELET: copy the LUT temperature into the solution ---*/ + Solution(iPoint,nDim+1) = FluidModel->GetTemperature(); + /*--- for FLAMELET: update the local temperature using LUT variables ---*/ + Temperature = Solution(iPoint,indices.Temperature()); + check_temp = SetTemperature(iPoint, Temperature); + /*--- Set the value of the density ---*/ const auto check_dens = SetDensity(iPoint, FluidModel->GetDensity()); @@ -84,7 +90,7 @@ bool CIncNSVariable::SetPrimVar(unsigned long iPoint, su2double eddy_visc, su2do /*--- Recompute the primitive variables ---*/ - Temperature = Solution(iPoint, nDim+1); + Temperature = Solution(iPoint, indices.Temperature()); SetTemperature(iPoint, Temperature); FluidModel->SetTDState_T(Temperature, scalar); SetDensity(iPoint, FluidModel->GetDensity()); diff --git a/SU2_CFD/src/variables/CSpeciesFlameletVariable.cpp b/SU2_CFD/src/variables/CSpeciesFlameletVariable.cpp new file mode 100644 index 00000000000..16536909bee --- /dev/null +++ b/SU2_CFD/src/variables/CSpeciesFlameletVariable.cpp @@ -0,0 +1,56 @@ +/*! + * \file CSpeciesFlameletVariable.cpp + * \brief Definition of the variable fields for the flamelet class. + * \author D. Mayer, T. Economon, N. Beishuizen + * \version 7.5.1 "Blackbird" + * + * SU2 Project Website: https://su2code.github.io + * + * The SU2 Project is maintained by the SU2 Foundation + * (http://su2foundation.org) + * + * Copyright 2012-2023, SU2 Contributors (cf. AUTHORS.md) + * + * SU2 is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public + * License as published by the Free Software Foundation; either + * version 2.1 of the License, or (at your option) any later version. + * + * SU2 is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with SU2. If not, see . + */ + +#include "../../include/variables/CSpeciesFlameletVariable.hpp" + +CSpeciesFlameletVariable::CSpeciesFlameletVariable(const su2double* species_inf, unsigned long npoint, + unsigned long ndim, unsigned long nvar, const CConfig* config) + : CSpeciesVariable(species_inf, npoint, ndim, nvar, config) { + for (unsigned long iPoint = 0; iPoint < nPoint; iPoint++) + for (unsigned long iVar = 0; iVar < nVar; iVar++) Solution(iPoint, iVar) = species_inf[iVar]; + + Solution_Old = Solution; + + /*--- Allocate and initialize solution for the dual time strategy ---*/ + bool dual_time = ((config->GetTime_Marching() == TIME_MARCHING::DT_STEPPING_1ST) || + (config->GetTime_Marching() == TIME_MARCHING::DT_STEPPING_2ND)); + + if (dual_time) { + Solution_time_n = Solution; + Solution_time_n1 = Solution; + } + + /*--- Allocate residual structures ---*/ + + Res_TruncError.resize(nPoint, nVar) = su2double(0.0); + + /* Allocate space for the source and scalars for visualization */ + + source_scalar.resize(nPoint, config->GetNScalars()) = su2double(0.0); + lookup_scalar.resize(nPoint, config->GetNLookups()) = su2double(0.0); + table_misses.resize(nPoint) = 0; +} diff --git a/TestCases/.gitignore b/TestCases/.gitignore index 315ff14f09e..a63f2d78d23 100644 --- a/TestCases/.gitignore +++ b/TestCases/.gitignore @@ -29,3 +29,6 @@ config_*.cfg *.tgz COPYING *.eqn + +# dragon files for flamelet models +*.drg diff --git a/TestCases/flamelet/01_laminar_premixed_ch4_flame_cfd/lam_prem_ch4_cfd.cfg b/TestCases/flamelet/01_laminar_premixed_ch4_flame_cfd/lam_prem_ch4_cfd.cfg new file mode 100644 index 00000000000..55ef51919de --- /dev/null +++ b/TestCases/flamelet/01_laminar_premixed_ch4_flame_cfd/lam_prem_ch4_cfd.cfg @@ -0,0 +1,135 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% SU2 configuration file % +% Case description: Laminar premixed flame stabilized on isothermal burner % +% Author: Nijso Beishuizen % +% Institution: Bosch Thermotechnology % +% Date: 08/09/2021 % +% File Version 7.5.1 "Blackbird" % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +SOLVER = INC_NAVIER_STOKES +KIND_TURB_MODEL= NONE +MATH_PROBLEM= DIRECT +RESTART_SOL = YES +% +% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% +% +INC_DENSITY_MODEL= VARIABLE +INC_ENERGY_EQUATION = YES +INC_DENSITY_INIT= 1.00 +INC_VELOCITY_INIT= (0.5, 0.0, 0.0 ) +INC_TEMPERATURE_INIT= 300.0 +INC_NONDIM= DIMENSIONAL +% +% -------------------- FLUID MODEL --------------------------------------- % +% +FLUID_MODEL= FLUID_FLAMELET +FILENAME_LUT= fgm_ch4.drg +% +% -------------------- SCALAR TRANSPORT ---------------------------------------% +% +KIND_SCALAR_MODEL= FLAMELET +DIFFUSIVITY_MODEL= FLAMELET +VISCOSITY_MODEL= FLAMELET +CONDUCTIVITY_MODEL= FLAMELET +FLAME_INIT= (0.0032, 0.00, 0.00, 1.0, 0.1, 0.0, 1.0e-3, 0.1) +% # progvar, enthalpy +SPECIES_INIT = (0.0, -210000, 0.0, 0.0) +CONV_NUM_METHOD_SPECIES= SCALAR_UPWIND +MUSCL_SPECIES= YES +SLOPE_LIMITER_SPECIES= NONE +TIME_DISCRE_SPECIES= EULER_IMPLICIT +% SCALAR CLIPPING +SPECIES_CLIPPING= YES +SPECIES_CLIPPING_MIN= 0 -2e6 0 0 +SPECIES_CLIPPING_MAX= 1 2e6 1 1 +% +MARKER_INLET_SPECIES = (inlet, 0, -193600.0, 0.0, 0.0) +CFL_REDUCTION_SPECIES= 1.0 +MARKER_SPECIES_STRONG_BC= (inlet, outlet, wall) +LOOKUP_NAMES = (MolarWeightMix, Conductivity, Cp, HeatRelease, Diffusivity) +USER_SCALAR_NAMES= (CO, NOx) +USER_SOURCE_NAMES = ( \ + ProdRatePos_Y-CO, ProdRateNegScaled_Y-CO ,\ + ProdRateTot_Y-NOx, NULL +) +% +% ---------------------- REFERENCE VALUE DEFINITION ---------------------------% +% +REF_ORIGIN_MOMENT_X = 0.25 +REF_ORIGIN_MOMENT_Y = 0.00 +REF_ORIGIN_MOMENT_Z = 0.00 +REF_LENGTH= 1.0 +REF_AREA= 1.0 +% +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +MARKER_ISOTHERMAL= (wall, 450) +MARKER_SYM= (symmetry, symmetry_2-part-fluid) +INC_INLET_TYPE= VELOCITY_INLET +MARKER_INLET = (inlet, 300.0, 0.20, 1.0, 0.0, 0.0) +INC_OUTLET_TYPE= PRESSURE_OUTLET +%INC_INLET_DAMPING = 0.05 +%INC_OUTLET_DAMPING = 0.05 +MARKER_OUTLET= (outlet, 0.0) +MARKER_PLOTTING= ( wall ) +MARKER_MONITORING= ( wall ) +MARKER_ANALYZE= ( inlet,outlet ) +MARKER_ANALYZE_AVERAGE = AREA +% +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +CFL_NUMBER= 350 +CFL_ADAPT= NO +ITER= 100 +OUTPUT_WRT_FREQ= 100 +% +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +LINEAR_SOLVER= FGMRES +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1E-03 +LINEAR_SOLVER_ITER= 10 +% +% -------------------------- MULTIGRID PARAMETERS -----------------------------% +% +MGLEVEL= 0 +% +% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% +% +CONV_NUM_METHOD_FLOW= FDS +MUSCL_FLOW= YES +SLOPE_LIMITER_FLOW = NONE +TIME_DISCRE_FLOW= EULER_IMPLICIT +% +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +CONV_RESIDUAL_MINVAL= -15 +CONV_STARTITER= 10 +CONV_CAUCHY_ELEMS= 100 +CONV_CAUCHY_EPS= 1E-6 +SCREEN_OUTPUT = INNER_ITER RMS_VELOCITY-X RMS_VELOCITY-Y RMS_PRESSURE RMS_PROGRESS_VARIABLE RMS_TOTAL_ENTHALPY RMS_CO RMS_NOx +HISTORY_OUTPUT = RMS_RES AERO_COEFF FLOW_COEFF FLOW_COEFF_SURF +VOLUME_OUTPUT = SOLUTION PRIMITIVE SOURCE RESIDUAL SENSITIVITY LOOKUP TIMESTEP +CONV_FIELD = RMS_PRESSURE +% +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +MESH_FORMAT= CGNS +MESH_FILENAME = mesh_structured.cgns +MESH_OUT_FILENAME= mesh_out.su2 +SOLUTION_FILENAME= solution +RESTART_FILENAME= restart +%OUTPUT_FILES = (RESTART,PARAVIEW,PARAVIEW_MULTIBLOCK) +OUTPUT_FILES = (RESTART) +TABULAR_FORMAT = CSV +CONV_FILENAME= history +VOLUME_FILENAME= ch4_flame_cfd +SURFACE_FILENAME= surface_flow +WRT_PERFORMANCE = NO +SCREEN_WRT_FREQ_INNER = 1 +SCREEN_WRT_FREQ_OUTER = 1 diff --git a/TestCases/flamelet/02_laminar_premixed_ch4_flame_hx_ad/lam_prem_ch4_hx_ad.cfg b/TestCases/flamelet/02_laminar_premixed_ch4_flame_hx_ad/lam_prem_ch4_hx_ad.cfg new file mode 100644 index 00000000000..884dd859268 --- /dev/null +++ b/TestCases/flamelet/02_laminar_premixed_ch4_flame_hx_ad/lam_prem_ch4_hx_ad.cfg @@ -0,0 +1,141 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Laminar premixed flame stabilized on isothermal burner % +% with cold wall (hx) downstream of the flame (adjoint run) % +% % +% Author a: Daniel Mayer % +% Institution a: Robert Bosch LLC % +% % +% Author b: Nijso Beishuizen % +% Institution b: Bosch Thermotechnology % +% % +% Date: 02/21/2023 % +% File Version 7.5.1 "Blackbird" % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +SOLVER= INC_NAVIER_STOKES +KIND_TURB_MODEL= NONE +% +OBJECTIVE_FUNCTION= CUSTOM_OBJFUNC +OBJECTIVE_WEIGHT= 1.0 +CUSTOM_OBJFUNC= avg_NOx +CUSTOM_OUTPUTS= 'avg_NOx : MassFlowAvg{SPECIES[3]}[outlet]' +% +% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% +% +INC_DENSITY_MODEL= VARIABLE +INC_DENSITY_INIT= 1.1766 +INC_VELOCITY_INIT= ( 0.2, 0, 0 ) +INC_ENERGY_EQUATION= YES +INC_TEMPERATURE_INIT= 300.0 +INC_NONDIM= DIMENSIONAL +% +% -------------------- FLUID PROPERTIES ------------------------------------- % +% +FLUID_MODEL= FLUID_FLAMELET +FILENAME_LUT= fgm_ch4.drg +CONDUCTIVITY_MODEL= FLAMELET +DIFFUSIVITY_MODEL= FLAMELET +KIND_SCALAR_MODEL= FLAMELET +VISCOSITY_MODEL= FLAMELET +% +FLAME_INIT= (0.0035, 0, 0, 3.0, 2.0, 0.0, 5.0e-4, 1.0) +% +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +MARKER_ISOTHERMAL= ( wall_burner, 450, wall_hx, 900 ) +MARKER_SYM= ( axis, symmetry_xm, symmetry_xp ) +% +INC_INLET_TYPE= VELOCITY_INLET +MARKER_INLET= (inlet, 300.0, 0.20, 1.0, 0.0, 0.0) +MARKER_INLET_SPECIES= (inlet, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) +MARKER_SPECIES_STRONG_BC= (inlet, wall_burner, wall_hx) + +% +INC_OUTLET_TYPE= PRESSURE_OUTLET +MARKER_OUTLET= ( outlet, 0.0 ) +% +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +CFL_NUMBER= 50 +CFL_REDUCTION_SPECIES= 1.0 +ITER= 11 +% +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +LINEAR_SOLVER= FGMRES +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1E-03 +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 ---------------------------------------% +% +USER_SCALAR_NAMES= (Y-CO, Y-NOx_tot) + +USER_SOURCE_NAMES= ( \ + ProdRatePos_Y-CO, ProdRateNegScaled_Y-CO ,\ + ProdRateTot_Y-NOx, zero) + +DIFFUSIVITY_CONSTANT= 0.00001 +% +CONV_NUM_METHOD_SPECIES= BOUNDED_SCALAR +MUSCL_SPECIES= NO +SLOPE_LIMITER_SPECIES= NONE +% +TIME_DISCRE_SPECIES= EULER_IMPLICIT +% +SPECIES_INIT= (0, 0, 0, 0, 0, 0, 0) +SPECIES_CLIPPING= NO +LOOKUP_NAMES= (MolarWeightMix, Conductivity, HeatRelease, Diffusivity, \ + ProdRatePos_Y-CO, ProdRateNegScaled_Y-CO ,\ + ProdRateTot_Y-NOx) + +% +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +CONV_FIELD= RMS_PRESSURE, RMS_VELOCITY-X, RMS_VELOCITY-Y +CONV_RESIDUAL_MINVAL= -18 +CONV_STARTITER= 10 +% +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +MESH_FILENAME= mesh_unstructured_hx.su2 +% +SCREEN_OUTPUT= INNER_ITER WALL_TIME RMS_ADJ_PRESSURE RMS_ADJ_VELOCITY-X RMS_ADJ_VELOCITY-Y RMS_ADJ_SPECIES_0 RMS_ADJ_SPECIES_1 RMS_ADJ_SPECIES_2 RMS_ADJ_SPECIES_3 +SCREEN_WRT_FREQ_INNER= 1 +% +HISTORY_OUTPUT= ITER RMS_RES LINSOL FLOW_COEFF_SURF CUSTOM COMBO +CONV_FILENAME= history +MARKER_ANALYZE= outlet +MARKER_ANALYZE_AVERAGE= MASSFLUX +% +%OUTPUT_FILES= RESTART, PARAVIEW, PARAVIEW_MULTIBLOCK +OUTPUT_FILES= RESTART, RESTART_ASCII +VOLUME_OUTPUT= RESIDUAL, PRIMITIVE +VOLUME_ADJ_FILENAME= ch4_flame_hx_ad + +OUTPUT_WRT_FREQ= 100 +% +GRAD_OBJFUNC_FILENAME= of_grad.csv +% +RESTART_SOL= YES +READ_BINARY_RESTART= NO +RESTART_FILENAME= restart +RESTART_ADJ_FILENAME= restart_adj +SOLUTION_FILENAME= solution +SOLUTION_ADJ_FILENAME= solution_adj +% +WRT_PERFORMANCE= YES +OUTPUT_PRECISION= 16 diff --git a/TestCases/flamelet/02_laminar_premixed_ch4_flame_hx_ad/lam_prem_ch4_hx_dot.cfg b/TestCases/flamelet/02_laminar_premixed_ch4_flame_hx_ad/lam_prem_ch4_hx_dot.cfg new file mode 100644 index 00000000000..14266a33b58 --- /dev/null +++ b/TestCases/flamelet/02_laminar_premixed_ch4_flame_hx_ad/lam_prem_ch4_hx_dot.cfg @@ -0,0 +1,251 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: Laminar premixed flame stabilized on isothermal burner % +% with cold wall (hx) downstream of the flame (adjoint run) % +% % +% Author a: Daniel Mayer % +% Institution a: Robert Bosch LLC % +% % +% Author b: Nijso Beishuizen % +% Institution b: Bosch Thermotechnology % +% % +% Date: 02/21/2023 % +% File Version 7.5.1 "Blackbird" % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +SOLVER= INC_NAVIER_STOKES +KIND_TURB_MODEL= NONE +% +OBJECTIVE_FUNCTION= CUSTOM_OBJFUNC +OBJECTIVE_WEIGHT= 1.0 +CUSTOM_OBJFUNC= avg_NOx +CUSTOM_OUTPUTS= 'avg_NOx : MassFlowAvg{SPECIES[3]}[outlet]' + +% +% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% +% +INC_DENSITY_MODEL= VARIABLE +INC_DENSITY_INIT= 1.1766 +INC_VELOCITY_INIT= ( 0.2, 0, 0 ) +INC_ENERGY_EQUATION= YES +INC_TEMPERATURE_INIT= 300.0 +INC_NONDIM= DIMENSIONAL +% +% -------------------- FLUID PROPERTIES ------------------------------------- % +% +FLUID_MODEL= FLUID_FLAMELET +FILENAME_LUT= fgm_ch4.drg +CONDUCTIVITY_MODEL= FLAMELET +DIFFUSIVITY_MODEL= FLAMELET +KIND_SCALAR_MODEL= FLAMELET +VISCOSITY_MODEL= FLAMELET +% +FLAME_INIT= (0.0035, 0, 0, 3, 2, 0, 5.0e-4, 1.0) +% +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +MARKER_ISOTHERMAL= ( wall_burner, 450, wall_hx, 900 ) +MARKER_SYM= ( axis, symmetry_xm, symmetry_xp, outlet ) +% +INC_INLET_TYPE= VELOCITY_INLET +MARKER_INLET= (inlet, 300.0, 0.20, 1.0, 0.0, 0.0) +MARKER_INLET_SPECIES= (inlet, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) +MARKER_SPECIES_STRONG_BC= (inlet, wall_burner, wall_hx) +% +% INC_OUTLET_TYPE= PRESSURE_OUTLET +% MARKER_OUTLET= ( outlet, 0.0 ) +% +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +CFL_NUMBER= 100 +CFL_REDUCTION_SPECIES= 1.0 +ITER= 11 +% +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +LINEAR_SOLVER= FGMRES +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1E-03 +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 ---------------------------------------% +% +USER_SCALAR_NAMES= (Y-CO, Y-NOx_tot) + +USER_SOURCE_NAMES= ( \ + ProdRatePos_Y-CO, ProdRateNegScaled_Y-CO ,\ + ProdRateTot_Y-NOx, zero) + +DIFFUSIVITY_CONSTANT= 0.00001 +% +CONV_NUM_METHOD_SPECIES= BOUNDED_SCALAR +MUSCL_SPECIES= NO +SLOPE_LIMITER_SPECIES= NONE +% +TIME_DISCRE_SPECIES= EULER_IMPLICIT +% +SPECIES_INIT= (0, 0, 0, 0, 0, 0, 0) +SPECIES_CLIPPING= NO +LOOKUP_NAMES= (MolarWeightMix, Conductivity, HeatRelease, Diffusivity, \ + ProdRatePos_Y-CO, ProdRateNegScaled_Y-CO ,\ + ProdRateTot_Y-NOx) + +% +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +CONV_FIELD= RMS_PRESSURE, RMS_VELOCITY-X, RMS_VELOCITY-Y +CONV_RESIDUAL_MINVAL= -18 +CONV_STARTITER= 10 +% +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +MESH_FILENAME= mesh_unstructured_hx.su2 +% +SCREEN_OUTPUT= INNER_ITER WALL_TIME RMS_ADJ_PRESSURE RMS_ADJ_VELOCITY-X RMS_ADJ_VELOCITY-Y RMS_ADJ_SPECIES_0 RMS_ADJ_SPECIES_1 RMS_ADJ_SPECIES_2 RMS_ADJ_SPECIES_3 +SCREEN_WRT_FREQ_INNER= 1 +% +HISTORY_OUTPUT= ITER RMS_RES LINSOL FLOW_COEFF_SURF +CONV_FILENAME= history +MARKER_ANALYZE= outlet +MARKER_ANALYZE_AVERAGE= MASSFLUX +% +%OUTPUT_FILES= RESTART, PARAVIEW, PARAVIEW_MULTIBLOCK +OUTPUT_FILES= RESTART +VOLUME_OUTPUT= RESIDUAL, PRIMITIVE +VOLUME_ADJ_FILENAME= ch4_flame_hx_ad + +OUTPUT_WRT_FREQ= 100 +% +GRAD_OBJFUNC_FILENAME= of_grad.csv +% +RESTART_SOL= YES +READ_BINARY_RESTART= NO +RESTART_FILENAME= restart +RESTART_ADJ_FILENAME= restart_adj +SOLUTION_FILENAME= solution +SOLUTION_ADJ_FILENAME= restart_adj +% +WRT_PERFORMANCE= YES +OUTPUT_PRECISION= 16 +% +% -------------------- FREE-FORM DEFORMATION PARAMETERS -----------------------% +% +% Tolerance of the Free-Form Deformation point inversion +FFD_TOLERANCE= 1E-10 +% +% Maximum number of iterations in the Free-Form Deformation point inversion +FFD_ITERATIONS= 500 +% +% FFD box definition: 3D case (FFD_BoxTag, X1, Y1, Z1, X2, Y2, Z2, X3, Y3, Z3, X4, Y4, Z4, +% X5, Y5, Z5, X6, Y6, Z6, X7, Y7, Z7, X8, Y8, Z8) +% 2D case (FFD_BoxTag, X1, Y1, 0.0, X2, Y2, 0.0, X3, Y3, 0.0, X4, Y4, 0.0, +% 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) +% Start at the lowest leftest corner and turn counter-clockwise +FFD_DEFINITION= (BOX, \ +0.065, 0.01, 0.0, \ +0.15, 0.01, 0.0 \ +0.15, 0.02, 0.0, \ +0.065, 0.02, 0.0, \ +0.0, 0.0 ,0.0, \ +0.0 ,0.0, 0.0, \ +0.0, 0.0, 0.0, \ +0.0, 0.0, 0.0 ) +% +% FFD box degree: 3D case (i_degree, j_degree, k_degree) +% 2D case (i_degree, j_degree, 0) +FFD_DEGREE= (6, 1, 0) +% +% Surface grid continuity at the intersection with the faces of the FFD boxes. +% To keep a particular level of surface continuity, SU2 automatically freezes the right +% number of control point planes (NO_DERIVATIVE, 1ST_DERIVATIVE, 2ND_DERIVATIVE, USER_INPUT) +FFD_CONTINUITY= USER_INPUT +% +% Definition of the FFD planes to be frozen in the FFD (x,y,z). +% Value from 0 FFD degree in that direction. Pick a value larger than degree if you do not want to fix any plane. +%FFD_FIX_I= (0,2,3) +%FFD_FIX_J= (0,2,3) +%FFD_FIX_K= (0,2,3) + +% ----------------------- DESIGN VARIABLE PARAMETERS --------------------------% +% + +DV_KIND= \ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D + +% +% Marker of the surface in which we are going apply the shape deformation +% NOTE: for deformation the outlet should be a MARKER_SYM to hinder the mesh being ripped apart. +DV_MARKER= ( wall_hx ) +% +% Parameters of the shape deformation +% - FFD_SETTING ( 1.0 ) +% - FFD_CONTROL_POINT ( FFD_BoxTag, i_Ind, j_Ind, k_Ind, x_Disp, y_Disp, z_Disp ) + +DV_PARAM= (BOX, 1, 0, 0, 1.0); \ +(BOX, 2, 0, 0, 1.0); \ +(BOX, 3, 0, 0, 1.0); \ +(BOX, 4, 0, 0, 1.0); \ +(BOX, 5, 0, 0, 1.0); \ +(BOX, 1, 1, 0, 1.0); \ +(BOX, 2, 1, 0, 1.0); \ +(BOX, 3, 1, 0, 1.0); \ +(BOX, 4, 1, 0, 1.0); \ +(BOX, 5, 1, 0, 1.0) + +DEFINITION_DV= ( 19, 1.0 | wall_hx | BOX, 0, 0, 1.0, 0.0 ) + +DV_VALUE= 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0 + +% +% ------------------------ GRID DEFORMATION PARAMETERS ------------------------% +% +% Linear solver or smoother for implicit formulations (FGMRES, RESTARTED_FGMRES, BCGSTAB) +DEFORM_LINEAR_SOLVER= FGMRES +% +% Preconditioner of the Krylov linear solver (ILU, LU_SGS, JACOBI) +DEFORM_LINEAR_SOLVER_PREC= ILU +% +% Number of smoothing iterations for mesh deformation +DEFORM_LINEAR_SOLVER_ITER= 1000 +% +% Number of nonlinear deformation iterations (surface deformation increments) +DEFORM_NONLINEAR_ITER= 1 +% +% Minimum residual criteria for the linear solver convergence of grid deformation +DEFORM_LINEAR_SOLVER_ERROR= 1E-14 +% +% Print the residuals during mesh deformation to the console (YES, NO) +DEFORM_CONSOLE_OUTPUT= YES +% +% Deformation coefficient (linear elasticity limits from -1.0 to 0.5, a larger +% value is also possible) +DEFORM_COEFF = 0.0 +% +% Type of element stiffness imposed for FEA mesh deformation (INVERSE_VOLUME, +% WALL_DISTANCE, CONSTANT_STIFFNESS) +DEFORM_STIFFNESS_TYPE= WALL_DISTANCE +% +% Deform the grid only close to the surface. It is possible to specify how much +% of the volumetric grid is going to be deformed in meters or inches (1E6 by default) +DEFORM_LIMIT = 1E6 diff --git a/TestCases/flamelet/03_laminar_premixed_ch4_flame_cht_cfd/lam_prem_ch4_cht_cfd_fluid.cfg b/TestCases/flamelet/03_laminar_premixed_ch4_flame_cht_cfd/lam_prem_ch4_cht_cfd_fluid.cfg new file mode 100644 index 00000000000..aff2d5561a8 --- /dev/null +++ b/TestCases/flamelet/03_laminar_premixed_ch4_flame_cht_cfd/lam_prem_ch4_cht_cfd_fluid.cfg @@ -0,0 +1,109 @@ + + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% + +SOLVER= INC_NAVIER_STOKES +KIND_TURB_MODEL= NONE + +MARKER_ANALYZE= (outlet) +MARKER_ANALYZE_AVERAGE= AREA + +MARKER_MONITORING= ( NONE ) +MARKER_PLOTTING= ( NONE ) + +WRT_ZONE_HIST= YES + +SCREEN_OUTPUT= INNER_ITER WALL_TIME RMS_PRESSURE RMS_PROGRESS_VARIABLE + +HISTORY_OUTPUT= RMS_RES AERO_COEFF FLOW_COEFF FLOW_COEFF_SURF + +VOLUME_OUTPUT= SOLUTION PRIMITIVE SOURCE RESIDUAL LOOKUP + +OUTPUT_FILES= (RESTART_ASCII) + +INNER_ITER= 1 + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% + +MARKER_SYM= ( symmetry_ym, symmetry_yp ) + +MARKER_ISOTHERMAL= ( wall_burner, 450.0 ) + +INLET_MATCHING_TOLERANCE= 1e-5 + +INC_INLET_TYPE= VELOCITY_INLET + +MARKER_INLET= (inlet, 300.0, 0.4, 1.0, 0.0, 0.0) +MARKER_INLET_SPECIES= (inlet, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) + +INC_OUTLET_TYPE= PRESSURE_OUTLET +MARKER_OUTLET= (outlet, 0.0) +MARKER_SPECIES_STRONG_BC= (inlet, outlet, wall_burner) + +% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% + +INC_DENSITY_MODEL= VARIABLE +INC_DENSITY_INIT= 1.1766 +INC_VELOCITY_INIT= ( 0.2, 0, 0 ) +INC_ENERGY_EQUATION= YES +INC_TEMPERATURE_INIT= 300.0 +INC_NONDIM= DIMENSIONAL + +% -------------------- FLUID MODEL --------------------------------------- % + +FLUID_MODEL= FLUID_FLAMELET +FILENAME_LUT= fgm_ch4.drg +FLAME_INIT= (0.004, 0.0, 0.0, 1.0, 0.0, 0.0, 0.2e-3, 1.0) + +% -------------------- SPECIES TRANSPORT ---------------------------------------% + +USER_SCALAR_NAMES= (Y-CO, Y-NOx_tot) + +USER_SOURCE_NAMES= ( \ + ProdRatePos_Y-CO, ProdRateNegScaled_Y-CO ,\ + ProdRateTot_Y-NOx, zero) + +DIFFUSIVITY_CONSTANT= 0.00001 + +CONV_NUM_METHOD_SPECIES= BOUNDED_SCALAR +MUSCL_SPECIES= NO +SLOPE_LIMITER_SPECIES= NONE + +TIME_DISCRE_SPECIES= EULER_IMPLICIT + +SPECIES_INIT= (0, 0, 0, 0, 0, 0, 0) +SPECIES_CLIPPING= NO +LOOKUP_NAMES= (MolarWeightMix, Conductivity, HeatRelease, Diffusivity, \ + ProdRatePos_Y-CO, ProdRateNegScaled_Y-CO ,\ + ProdRateTot_Y-NOx) + +% --------------------------- VISCOSITY MODEL ---------------------------------% + +CONDUCTIVITY_MODEL= FLAMELET +DIFFUSIVITY_MODEL= FLAMELET +KIND_SCALAR_MODEL= FLAMELET +VISCOSITY_MODEL= FLAMELET + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% + +%NUM_METHOD_GRAD= GREEN_GAUSS +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES + +CFL_NUMBER= 25 +CFL_ADAPT= NO + +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% + +LINEAR_SOLVER= BCGSTAB +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1E-10 +LINEAR_SOLVER_ITER= 20 + +% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% + +CONV_NUM_METHOD_FLOW= FDS +MUSCL_FLOW= NO +SLOPE_LIMITER_FLOW= NONE +VENKAT_LIMITER_COEFF= 10.0 +JST_SENSOR_COEFF= ( 0.5, 0.02 ) +TIME_DISCRE_FLOW= EULER_IMPLICIT diff --git a/TestCases/flamelet/03_laminar_premixed_ch4_flame_cht_cfd/lam_prem_ch4_cht_cfd_master.cfg b/TestCases/flamelet/03_laminar_premixed_ch4_flame_cht_cfd/lam_prem_ch4_cht_cfd_master.cfg new file mode 100644 index 00000000000..1ca70ab9d79 --- /dev/null +++ b/TestCases/flamelet/03_laminar_premixed_ch4_flame_cht_cfd/lam_prem_ch4_cht_cfd_master.cfg @@ -0,0 +1,55 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: flow through a heat exchanger, CHT setup % +% Author: Nijso Beishuizen, Ole Burghardt % +% Institution: Bosch Thermotechniek BV, Technical University of Kaiserslautern % +% Date: 24/11/2020 % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +SOLVER= MULTIPHYSICS + +CONFIG_LIST= ( lam_prem_ch4_cht_cfd_fluid.cfg, lam_prem_ch4_cht_cfd_solid.cfg ) + +%%%%%%%%%%%%%%%%%% +% Mesh stuff % +%%%%%%%%%%%%%%%%%% + +MESH_FORMAT= SU2 +MULTIZONE_MESH= YES + +% Restart flow input file +SOLUTION_FILENAME= solution + +% define the interface between the zones +MARKER_ZONE_INTERFACE= ( cht_interface_fluid, cht_interface_solid ) +MARKER_CHT_INTERFACE= ( cht_interface_fluid, cht_interface_solid ) + +CHT_COUPLING_METHOD= DIRECT_TEMPERATURE_ROBIN_HEATFLUX + +READ_BINARY_RESTART= NO +RESTART_SOL= YES + +RESTART_FILENAME= restart +RESTART_ADJ_FILENAME= restart_adj + +TIME_DOMAIN= NO +OUTER_ITER= 10 + +OUTPUT_WRT_FREQ= 500 +SCREEN_WRT_FREQ_INNER= 1 +SCREEN_WRT_FREQ_OUTER= 1 + +WRT_VOLUME_OVERWRITE= YES + +WRT_ZONE_CONV= YES +CONV_RESIDUAL_MINVAL= -20 + +VOLUME_FILENAME= fluid + +CONV_FILENAME= history + +% --------------------------- Optimization Parameters --------------------------% + +MESH_FILENAME= mesh_unstructured_cht.su2 diff --git a/TestCases/flamelet/03_laminar_premixed_ch4_flame_cht_cfd/lam_prem_ch4_cht_cfd_solid.cfg b/TestCases/flamelet/03_laminar_premixed_ch4_flame_cht_cfd/lam_prem_ch4_cht_cfd_solid.cfg new file mode 100644 index 00000000000..197980c4899 --- /dev/null +++ b/TestCases/flamelet/03_laminar_premixed_ch4_flame_cht_cfd/lam_prem_ch4_cht_cfd_solid.cfg @@ -0,0 +1,57 @@ + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% + +SOLVER= HEAT_EQUATION + +MARKER_PLOTTING= ( NONE ) +MARKER_MONITORING= ( NONE ) + +WRT_ZONE_HIST= YES +OUTPUT_FILES= ( RESTART_ASCII ) + +INC_NONDIM= DIMENSIONAL + +INNER_ITER= 1 + +% ---------------- (SOLIDS) CONDUCTION CONDITION DEFINITION -------------------% + +FREESTREAM_TEMPERATURE= 450.0 + +MATERIAL_DENSITY= 2700 + +THERMAL_CONDUCTIVITY_CONSTANT=200 + +SPECIFIC_HEAT_CP= 900 + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% + +MARKER_ISOTHERMAL= ( wall_water, 350.0 ) + +% -------------------- HEAT NUMERICAL METHOD DEFINITION -----------------------% + +TIME_DISCRE_HEAT= EULER_IMPLICIT + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% + +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +CFL_NUMBER= 250 +CFL_ADAPT= NO + +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% + +LINEAR_SOLVER= BCGSTAB +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1e-16 +LINEAR_SOLVER_ITER= 100 + +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% + +TABULAR_FORMAT= CSV +CONV_FILENAME= history +VOLUME_FILENAME= solid +VOLUME_ADJ_FILENAME= adjoint +VALUE_OBJFUNC_FILENAME= of_eval.dat +SURFACE_FILENAME= surface_solid +SURFACE_ADJ_FILENAME= surface_adjoint + + diff --git a/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_ad_fluid.cfg b/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_ad_fluid.cfg new file mode 100644 index 00000000000..b100aa3fcb3 --- /dev/null +++ b/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_ad_fluid.cfg @@ -0,0 +1,110 @@ + + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% + +SOLVER= INC_NAVIER_STOKES +KIND_TURB_MODEL= NONE + +MARKER_ANALYZE= (outlet) +MARKER_ANALYZE_AVERAGE= AREA + +MARKER_MONITORING= ( NONE ) +MARKER_PLOTTING= ( NONE ) + +WRT_ZONE_HIST= YES + +%SCREEN_OUTPUT= INNER_ITER WALL_TIME RMS_PRESSURE RMS_VELOCITY-X RMS_VELOCITY-Y RMS_TKE RMS_DISSIPATION LINSOL_ITER LINSOL_RESIDUAL RMS_ADJ_PRESSURE RMS_ADJ_VELOCITY-X RMS_ADJ_VELOCITY-Y + +SCREEN_OUTPUT= INNER_ITER WALL_TIME RMS_ADJ_PRESSURE RMS_ADJ_VELOCITY-X RMS_ADJ_VELOCITY-Y RMS_ADJ_PROGRESS_VARIABLE RMS_ADJ_TOTAL_ENTHALPY RMS_ADJ_Y-CO RMS_ADJ_Y-NOx_tot + +HISTORY_OUTPUT= RMS_RES AERO_COEFF FLOW_COEFF FLOW_COEFF_SURF + +VOLUME_OUTPUT= SOLUTION PRIMITIVE SOURCE RESIDUAL LOOKUP + +OUTPUT_FILES= ( RESTART_ASCII) + +INNER_ITER= 1 + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% + +MARKER_SYM= ( symmetry_ym, symmetry_yp ) + +MARKER_ISOTHERMAL= ( wall_burner, 450.0 ) + +INLET_MATCHING_TOLERANCE= 1e-5 + +INC_INLET_TYPE= VELOCITY_INLET + +MARKER_INLET= (inlet, 300.0, 0.4, 1.0, 0.0, 0.0) +MARKER_INLET_SPECIES= (inlet, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) + +INC_OUTLET_TYPE= PRESSURE_OUTLET +MARKER_OUTLET= (outlet, 0.0) +MARKER_SPECIES_STRONG_BC= (inlet, wall_burner) + +% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% + +INC_DENSITY_MODEL= VARIABLE +INC_DENSITY_INIT= 1.1766 +INC_VELOCITY_INIT= ( 0.2, 0, 0 ) +INC_ENERGY_EQUATION= YES +INC_TEMPERATURE_INIT= 300.0 +INC_NONDIM= DIMENSIONAL + +% -------------------- FLUID MODEL --------------------------------------- % + +FLUID_MODEL= FLUID_FLAMELET +FILENAME_LUT= fgm_ch4.drg +FLAME_INIT= (0.004, 0.0, 0.0, 1, 0, 0, 0.2e-3, 1) + +% -------------------- SPECIES TRANSPORT ---------------------------------------% + +USER_SCALAR_NAMES= (Y-CO, Y-NOx_tot) + +USER_SOURCE_NAMES= ( \ + ProdRatePos_Y-CO, ProdRateNegScaled_Y-CO ,\ + ProdRateTot_Y-NOx, zero) + +DIFFUSIVITY_CONSTANT= 0.00001 + +CONV_NUM_METHOD_SPECIES= BOUNDED_SCALAR +MUSCL_SPECIES= NO +SLOPE_LIMITER_SPECIES= NONE + +TIME_DISCRE_SPECIES= EULER_IMPLICIT + +SPECIES_INIT= (0, 0, 0, 0, 0, 0, 0) +SPECIES_CLIPPING= NO +LOOKUP_NAMES= (MolarWeightMix, Conductivity, HeatRelease, Diffusivity, \ + ProdRatePos_Y-CO, ProdRateNegScaled_Y-CO ,\ + ProdRateTot_Y-NOx) + +% --------------------------- VISCOSITY MODEL ---------------------------------% + +CONDUCTIVITY_MODEL= FLAMELET +DIFFUSIVITY_MODEL= FLAMELET +KIND_SCALAR_MODEL= FLAMELET +VISCOSITY_MODEL= FLAMELET + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% + +%NUM_METHOD_GRAD= GREEN_GAUSS +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +CFL_NUMBER= 25 +CFL_ADAPT= NO + +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% + +LINEAR_SOLVER= BCGSTAB +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1E-10 +LINEAR_SOLVER_ITER= 20 + +% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% + +CONV_NUM_METHOD_FLOW= FDS +MUSCL_FLOW= NO +SLOPE_LIMITER_FLOW= NONE +VENKAT_LIMITER_COEFF= 10.0 +JST_SENSOR_COEFF= ( 0.5, 0.02 ) +TIME_DISCRE_FLOW= EULER_IMPLICIT diff --git a/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_ad_master.cfg b/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_ad_master.cfg new file mode 100644 index 00000000000..e3e9a0597c1 --- /dev/null +++ b/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_ad_master.cfg @@ -0,0 +1,97 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: flow through a heat exchanger, CHT setup % +% Author: Nijso Beishuizen, Ole Burghardt % +% Institution: Bosch Thermotechniek BV, Technical University of Kaiserslautern % +% Date: 24/11/2020 % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +SOLVER= MULTIPHYSICS + +CONFIG_LIST= ( lam_prem_ch4_cht_ad_fluid.cfg, lam_prem_ch4_cht_ad_solid.cfg ) + + +%%%%%%%%%%%%%%%%%% +% Mesh stuff % +%%%%%%%%%%%%%%%%%% + + +MESH_FILENAME= mesh_unstructured_cht.su2 +% MESH_OUT_FILENAME= mesh_def.su2 + +MESH_FORMAT= SU2 +MULTIZONE_MESH= YES + +% Restart flow input file +SOLUTION_FILENAME= solution +READ_BINARY_RESTART= NO +RESTART_SOL= YES +RESTART_FILENAME= restart +RESTART_ADJ_FILENAME= restart_adj +SOLUTION_ADJ_FILENAME= solution_adj + +% define the interface between the zones +MARKER_ZONE_INTERFACE= ( cht_interface_fluid, cht_interface_solid ) +MARKER_CHT_INTERFACE= ( cht_interface_fluid, cht_interface_solid ) + +CHT_COUPLING_METHOD= DIRECT_TEMPERATURE_ROBIN_HEATFLUX + +TIME_DOMAIN= NO +OUTER_ITER = 1100 + +OUTPUT_WRT_FREQ= 50 +SCREEN_WRT_FREQ_INNER= 1 +SCREEN_WRT_FREQ_OUTER= 1 + +WRT_VOLUME_OVERWRITE= YES + +WRT_ZONE_CONV= YES +CONV_RESIDUAL_MINVAL= -20 + +VOLUME_FILENAME= fluid + +CONV_FILENAME= history + +% --------------------------- Optimization Parameters --------------------------% + +GRAD_OBJFUNC_FILENAME= of_grad.csv + +OBJECTIVE_FUNCTION= SURFACE_STATIC_TEMPERATURE +OBJECTIVE_WEIGHT= 1.0 + +DV_KIND= \ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D + +DV_MARKER= ( cht_interface_fluid, cht_interface_solid ) + +DV_PARAM= \ + (BOX, 0, 0, 1.0, 1.0); \ + (BOX, 2, 0, 1.0, 1.0); \ + (BOX, 2, 2, 1.0, 1.0); \ + (BOX, 0, 2, 1.0, 1.0); \ +% + (BOX, 1, 0, 0.0, 1.0); \ + (BOX, 2, 1, 1.0, 0.0); \ + (BOX, 1, 2, 0.0, 1.0); \ + (BOX, 0, 1, 1.0, 0.0) + +FFD_DEFINITION= (BOX, 0.0089, 0.0006, 0.0, 0.018, 0.0006, 0.0, 0.018, 0.0018, 0.0, 0.0089, 0.0018, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) +FFD_DEGREE= ( 2, 1, 0) +FFD_CONTINUITY= 1ST_DERIVATIVE +FFD_FIX_I= (100) +FFD_FIX_J= (100) +FFD_ITERATIONS= 1000 +FFD_TOLERANCE= 1E-10 + +DEFORM_STIFFNESS_TYPE= INVERSE_VOLUME + +DV_VALUE=0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 diff --git a/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_ad_solid.cfg b/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_ad_solid.cfg new file mode 100644 index 00000000000..197980c4899 --- /dev/null +++ b/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_ad_solid.cfg @@ -0,0 +1,57 @@ + +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% + +SOLVER= HEAT_EQUATION + +MARKER_PLOTTING= ( NONE ) +MARKER_MONITORING= ( NONE ) + +WRT_ZONE_HIST= YES +OUTPUT_FILES= ( RESTART_ASCII ) + +INC_NONDIM= DIMENSIONAL + +INNER_ITER= 1 + +% ---------------- (SOLIDS) CONDUCTION CONDITION DEFINITION -------------------% + +FREESTREAM_TEMPERATURE= 450.0 + +MATERIAL_DENSITY= 2700 + +THERMAL_CONDUCTIVITY_CONSTANT=200 + +SPECIFIC_HEAT_CP= 900 + +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% + +MARKER_ISOTHERMAL= ( wall_water, 350.0 ) + +% -------------------- HEAT NUMERICAL METHOD DEFINITION -----------------------% + +TIME_DISCRE_HEAT= EULER_IMPLICIT + +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% + +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +CFL_NUMBER= 250 +CFL_ADAPT= NO + +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% + +LINEAR_SOLVER= BCGSTAB +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1e-16 +LINEAR_SOLVER_ITER= 100 + +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% + +TABULAR_FORMAT= CSV +CONV_FILENAME= history +VOLUME_FILENAME= solid +VOLUME_ADJ_FILENAME= adjoint +VALUE_OBJFUNC_FILENAME= of_eval.dat +SURFACE_FILENAME= surface_solid +SURFACE_ADJ_FILENAME= surface_adjoint + + diff --git a/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_dot_master.cfg b/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_dot_master.cfg new file mode 100644 index 00000000000..4a3e831b6c6 --- /dev/null +++ b/TestCases/flamelet/04_laminar_premixed_ch4_flame_cht_ad/lam_prem_ch4_cht_dot_master.cfg @@ -0,0 +1,97 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% % +% SU2 configuration file % +% Case description: flow through a heat exchanger, CHT setup % +% Author: Nijso Beishuizen, Ole Burghardt % +% Institution: Bosch Thermotechniek BV, Technical University of Kaiserslautern % +% Date: 24/11/2020 % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +SOLVER= MULTIPHYSICS + +CONFIG_LIST= ( lam_prem_ch4_cht_ad_fluid.cfg, lam_prem_ch4_cht_ad_solid.cfg ) + + +%%%%%%%%%%%%%%%%%% +% Mesh stuff % +%%%%%%%%%%%%%%%%%% + + +MESH_FILENAME= mesh_unstructured_cht.su2 +% MESH_OUT_FILENAME= mesh_def.su2 + +MESH_FORMAT= SU2 +MULTIZONE_MESH= YES + +% Restart flow input file +SOLUTION_FILENAME= solution +READ_BINARY_RESTART= NO +RESTART_SOL= YES +RESTART_FILENAME= restart +RESTART_ADJ_FILENAME= restart_adj +SOLUTION_ADJ_FILENAME= restart_adj + +% define the interface between the zones +MARKER_ZONE_INTERFACE= ( cht_interface_fluid, cht_interface_solid ) +MARKER_CHT_INTERFACE= ( cht_interface_fluid, cht_interface_solid ) + +CHT_COUPLING_METHOD= DIRECT_TEMPERATURE_ROBIN_HEATFLUX + +TIME_DOMAIN= NO +OUTER_ITER = 11 + +OUTPUT_WRT_FREQ= 500 +SCREEN_WRT_FREQ_INNER= 1 +SCREEN_WRT_FREQ_OUTER= 1 + +WRT_VOLUME_OVERWRITE= YES + +WRT_ZONE_CONV= YES +CONV_RESIDUAL_MINVAL= -20 + +VOLUME_FILENAME= fluid + +CONV_FILENAME= history + +% --------------------------- Optimization Parameters --------------------------% + +GRAD_OBJFUNC_FILENAME= of_grad.csv + +OBJECTIVE_FUNCTION= SURFACE_STATIC_TEMPERATURE +OBJECTIVE_WEIGHT= 1.0 + +DV_KIND= \ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D,\ +FFD_CONTROL_POINT_2D + +DV_MARKER= ( cht_interface_fluid, cht_interface_solid ) + +DV_PARAM= \ + (BOX, 0, 0, 1.0, 1.0); \ + (BOX, 2, 0, 1.0, 1.0); \ + (BOX, 2, 2, 1.0, 1.0); \ + (BOX, 0, 2, 1.0, 1.0); \ +% + (BOX, 1, 0, 0.0, 1.0); \ + (BOX, 2, 1, 1.0, 0.0); \ + (BOX, 1, 2, 0.0, 1.0); \ + (BOX, 0, 1, 1.0, 0.0) + +FFD_DEFINITION= (BOX, 0.0089, 0.0006, 0.0, 0.018, 0.0006, 0.0, 0.018, 0.0018, 0.0, 0.0089, 0.0018, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0) +FFD_DEGREE= ( 2, 1, 0) +FFD_CONTINUITY= 1ST_DERIVATIVE +FFD_FIX_I= (100) +FFD_FIX_J= (100) +FFD_ITERATIONS= 1000 +FFD_TOLERANCE= 1E-10 + +DEFORM_STIFFNESS_TYPE= INVERSE_VOLUME + +DV_VALUE=0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0 diff --git a/TestCases/flamelet/05_laminar_premixed_ch4_flame_cfd_axi/lam_prem_ch4_cfd_axi.cfg b/TestCases/flamelet/05_laminar_premixed_ch4_flame_cfd_axi/lam_prem_ch4_cfd_axi.cfg new file mode 100644 index 00000000000..6fa6662e253 --- /dev/null +++ b/TestCases/flamelet/05_laminar_premixed_ch4_flame_cfd_axi/lam_prem_ch4_cfd_axi.cfg @@ -0,0 +1,136 @@ +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% SU2 configuration file % +% Case description: Laminar premixed flame stabilized on isothermal burner % +% axisymmetric setup +% Author: Nijso Beishuizen % +% Institution: Bosch Thermotechnology % +% Date: 08/09/2021 % +% File Version 7.5.0 "Falcon", branch: feature_new_flamelet % +% % +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% +% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------% +% +SOLVER= INC_NAVIER_STOKES +KIND_TURB_MODEL= NONE +MATH_PROBLEM= DIRECT +RESTART_SOL= YES +% +% ---------------- INCOMPRESSIBLE FLOW CONDITION DEFINITION -------------------% +% +INC_DENSITY_MODEL= VARIABLE +INC_ENERGY_EQUATION= YES +INC_DENSITY_INIT= 1.00 +INC_VELOCITY_INIT= (0.5, 0.0, 0.0 ) +INC_TEMPERATURE_INIT= 300.0 +INC_NONDIM= DIMENSIONAL +% +% +AXISYMMETRIC= YES +% -------------------- FLUID MODEL --------------------------------------- % +% +FLUID_MODEL= FLUID_FLAMELET +FILENAME_LUT= fgm_ch4.drg +% +% -------------------- SCALAR TRANSPORT ---------------------------------------% +% +KIND_SCALAR_MODEL= FLAMELET +DIFFUSIVITY_MODEL= FLAMELET +VISCOSITY_MODEL= FLAMELET +CONDUCTIVITY_MODEL= FLAMELET +FLAME_INIT= (0.0032, 0.0, 0.0, 1, 0.1, 0, 1e-3, 0.1) +% # progvar, enthalpy +SPECIES_INIT= (0.0, -210000, 0.0, 0.0) +CONV_NUM_METHOD_SPECIES= BOUNDED_SCALAR +%CONV_NUM_METHOD_SPECIES= SCALAR_UPWIND +MUSCL_SPECIES= YES +SLOPE_LIMITER_SPECIES= NONE +TIME_DISCRE_SPECIES= EULER_IMPLICIT +% SCALAR CLIPPING +SPECIES_CLIPPING= YES +SPECIES_CLIPPING_MIN= 0 -2e6 0 0 +SPECIES_CLIPPING_MAX= 1 2e6 1 1 +% +MARKER_INLET_SPECIES= (inlet, 0, -193600.0, 0.0, 0.0) +CFL_REDUCTION_SPECIES= 1.0 +MARKER_SPECIES_STRONG_BC= (inlet, outlet, wall) +LOOKUP_NAMES= (MolarWeightMix, Conductivity, Cp, HeatRelease, Diffusivity) +USER_SCALAR_NAMES= (CO, NOx) +USER_SOURCE_NAMES= ( \ + ProdRatePos_Y-CO, ProdRateNegScaled_Y-CO ,\ + ProdRateTot_Y-NOx, NULL ) +% +% ---------------------- REFERENCE VALUE DEFINITION ---------------------------% +% +REF_ORIGIN_MOMENT_X= 0.25 +REF_ORIGIN_MOMENT_Y= 0.00 +REF_ORIGIN_MOMENT_Z= 0.00 +REF_LENGTH= 1.0 +REF_AREA= 1.0 +% +% -------------------- BOUNDARY CONDITION DEFINITION --------------------------% +% +MARKER_ISOTHERMAL= (wall, 450) +MARKER_SYM= (symmetry, symmetry_2-part-fluid) +INC_INLET_TYPE= VELOCITY_INLET +MARKER_INLET= (inlet, 300.0, 0.20, 1.0, 0.0, 0.0) +INC_OUTLET_TYPE= PRESSURE_OUTLET +MARKER_OUTLET= (outlet, 0.0) +MARKER_PLOTTING= ( wall ) +MARKER_MONITORING= ( wall ) +MARKER_ANALYZE= ( inlet,outlet ) +MARKER_ANALYZE_AVERAGE= AREA +% +% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------% +% +NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES +CFL_NUMBER= 20 +CFL_ADAPT= NO +ITER= 10000 +OUTPUT_WRT_FREQ= 100 +% +% ------------------------ LINEAR SOLVER DEFINITION ---------------------------% +% +LINEAR_SOLVER= FGMRES +LINEAR_SOLVER_PREC= ILU +LINEAR_SOLVER_ERROR= 1E-10 +LINEAR_SOLVER_ITER= 25 +% +% -------------------------- MULTIGRID PARAMETERS -----------------------------% +% +MGLEVEL= 0 +% +% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------% +% +CONV_NUM_METHOD_FLOW= FDS +MUSCL_FLOW= YES +SLOPE_LIMITER_FLOW= NONE +TIME_DISCRE_FLOW= EULER_IMPLICIT +% +% --------------------------- CONVERGENCE PARAMETERS --------------------------% +% +CONV_RESIDUAL_MINVAL= -15 +CONV_STARTITER= 10 +CONV_CAUCHY_ELEMS= 100 +CONV_CAUCHY_EPS= 1E-6 +SCREEN_OUTPUT= INNER_ITER WALL_TIME RMS_VELOCITY-X RMS_VELOCITY-Y RMS_PRESSURE RMS_PROGRESS_VARIABLE RMS_TOTAL_ENTHALPY +HISTORY_OUTPUT= RMS_RES AERO_COEFF FLOW_COEFF FLOW_COEFF_SURF +VOLUME_OUTPUT= SOLUTION PRIMITIVE SOURCE RESIDUAL SENSITIVITY LOOKUP TIMESTEP +CONV_FIELD= RMS_PRESSURE +% +% ------------------------- INPUT/OUTPUT INFORMATION --------------------------% +% +MESH_FORMAT= SU2 +MESH_FILENAME= mesh_axi.su2 +MESH_OUT_FILENAME= mesh_out.su2 +SOLUTION_FILENAME= solution +RESTART_FILENAME= restart +%OUTPUT_FILES= (RESTART,PARAVIEW,PARAVIEW_MULTIBLOCK) +OUTPUT_FILES= (RESTART) +TABULAR_FORMAT= CSV +CONV_FILENAME= history +VOLUME_FILENAME= flow +SURFACE_FILENAME= surface_flow +WRT_PERFORMANCE= NO +SCREEN_WRT_FREQ_INNER= 1 +SCREEN_WRT_FREQ_OUTER= 1 diff --git a/TestCases/flamelet/README.md b/TestCases/flamelet/README.md new file mode 100644 index 00000000000..8608d466f0d --- /dev/null +++ b/TestCases/flamelet/README.md @@ -0,0 +1,8 @@ +# Flamelet regression cases + +## laminar_premixed_flame +1. A laminar premixed flame stabilized on an isothermal burner with a fixed wall temperature using a low resolution look-up table for all thermo-chemical quantities. +2. A laminar premixed flame including a heat exchanger, AD test. +3. A laminar premixed flame including a heat exchanger (Conjugate Heat Transfer). +4. A laminar premixed flame including a heat exchanger (Conjugate Heat Transfer), AD test. +5. An axisymmetric laminar premixed flame. diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index 95c486de561..50dd5b013ce 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -38,6 +38,28 @@ def main(): test_list = [] + ####################### + ### Flamelet solver ### + ####################### + + # 2D planar laminar premixed methane flame on isothermal burner (restart) + cfd_flamelet_ch4 = TestCase('cfd_flamelet_ch4') + cfd_flamelet_ch4.cfg_dir = "flamelet/01_laminar_premixed_ch4_flame_cfd" + cfd_flamelet_ch4.cfg_file = "lam_prem_ch4_cfd.cfg" + cfd_flamelet_ch4.test_iter = 10 + cfd_flamelet_ch4.test_vals = [-15.313265, -15.180884, -15.291808, -8.488238, -15.010141, -15.920950] + cfd_flamelet_ch4.new_output = True + test_list.append(cfd_flamelet_ch4) + + # axisymmetric 2D planar laminar premixed methane flame on isothermal burner (restart) + cfd_flamelet_ch4_axi = TestCase('cfd_flamelet_ch4_axi') + cfd_flamelet_ch4_axi.cfg_dir = "flamelet/05_laminar_premixed_ch4_flame_cfd_axi" + cfd_flamelet_ch4_axi.cfg_file = "lam_prem_ch4_cfd_axi.cfg" + cfd_flamelet_ch4_axi.test_iter = 10 + cfd_flamelet_ch4_axi.test_vals = [-11.054149, -12.276393, -11.299388, -13.877670, -6.291548] + cfd_flamelet_ch4_axi.new_output = True + test_list.append(cfd_flamelet_ch4_axi) + ######################### ## NEMO solver ### ######################### @@ -1502,6 +1524,25 @@ def main(): pass_list = [ test.run_test() for test in test_list ] + ###################################### + ### RUN CHT TEST WITH FILEDIFF ### + ###################################### + + # 2D planar laminar premixed methane flame on isothermal burner with conjugate heat transfer in cooling fin (restart) + cfd_flamelet_ch4_cht = TestCase('cfd_flamelet_ch4_cht') + cfd_flamelet_ch4_cht.cfg_dir = "flamelet/03_laminar_premixed_ch4_flame_cht_cfd" + cfd_flamelet_ch4_cht.cfg_file = "lam_prem_ch4_cht_cfd_master.cfg" + cfd_flamelet_ch4_cht.test_iter = 10 + cfd_flamelet_ch4_cht.command = TestCase.Command("mpirun -n 2", "SU2_CFD") + cfd_flamelet_ch4_cht.timeout = 1600 + cfd_flamelet_ch4_cht.reference_file = "restart_0.csv.ref" + cfd_flamelet_ch4_cht.test_file = "restart_0.csv" + cfd_flamelet_ch4_cht.multizone = True + cfd_flamelet_ch4_cht.comp_threshold = 1e-6 + cfd_flamelet_ch4_cht.tol_file_percent = 0.1 + pass_list.append(cfd_flamelet_ch4_cht.run_filediff()) + test_list.append(cfd_flamelet_ch4_cht) + ###################################### ### RUN SU2_SOL TESTS ### ###################################### diff --git a/TestCases/parallel_regression_AD.py b/TestCases/parallel_regression_AD.py index 436715d9866..d7415431119 100644 --- a/TestCases/parallel_regression_AD.py +++ b/TestCases/parallel_regression_AD.py @@ -340,6 +340,77 @@ def main(): pass_list = [ test.run_test() for test in test_list ] + ################################## + ### Disc. adj. flamelet solver ### + ################################## + + # 2D planar laminar premixed flame on isothermal burner (restart) + discadj_flamelet_ch4_hx = TestCase('discadj_flamelet_ch4_hx') + discadj_flamelet_ch4_hx.command = TestCase.Command("mpirun -n 2", "SU2_CFD_AD") + discadj_flamelet_ch4_hx.cfg_dir = "flamelet/02_laminar_premixed_ch4_flame_hx_ad" + discadj_flamelet_ch4_hx.cfg_file = "lam_prem_ch4_hx_ad.cfg" + discadj_flamelet_ch4_hx.multizone = True + discadj_flamelet_ch4_hx.test_iter = 10 + discadj_flamelet_ch4_hx.timeout = 20000 + discadj_flamelet_ch4_hx.reference_file = "restart_adj_custom.csv.ref" + discadj_flamelet_ch4_hx.test_file = "restart_adj_custom.csv" + discadj_flamelet_ch4_hx.comp_threshold = 1e-6 + discadj_flamelet_ch4_hx.tol_file_percent = 0.1 + pass_list.append(discadj_flamelet_ch4_hx.run_filediff()) + test_list.append(discadj_flamelet_ch4_hx) + + # 2D planar laminar premixed flame on isothermal burner with conjugate heat transfer (restart) + discadj_flamelet_ch4_cht = TestCase('discadj_flamelet_ch4_cht') + discadj_flamelet_ch4_cht.command = TestCase.Command("mpirun -n 2", "SU2_CFD_AD") + discadj_flamelet_ch4_cht.cfg_dir = "flamelet/04_laminar_premixed_ch4_flame_cht_ad" + discadj_flamelet_ch4_cht.cfg_file = "lam_prem_ch4_cht_ad_master.cfg" + discadj_flamelet_ch4_cht.multizone = True + discadj_flamelet_ch4_cht.test_iter = 5 + discadj_flamelet_ch4_cht.reference_file = "restart_adj_T_0.csv.ref" + discadj_flamelet_ch4_cht.test_file = "restart_adj_T_0.csv" + discadj_flamelet_ch4_cht.comp_threshold = 1e-6 + discadj_flamelet_ch4_cht.tol_file_percent = 0.1 + discadj_flamelet_ch4_cht.timeout = 20000 + pass_list.append(discadj_flamelet_ch4_cht.run_filediff()) + test_list.append(discadj_flamelet_ch4_cht) + + ################################################ + ### Gradient check (dot) for flamelet solver ### + ################################################ + + # 2D planar laminar premixed flame on isothermal burner (restart) + # This test restarts on the output of test discadj_flamelet_ch4_hx and + # will only pass if test discadj_flamelet_ch4_hx passes. + dot_flamelet_ch4_hx = TestCase('dot_flamelet_ch4_hx') + dot_flamelet_ch4_hx.cfg_dir = "flamelet/02_laminar_premixed_ch4_flame_hx_ad" + dot_flamelet_ch4_hx.cfg_file = "lam_prem_ch4_hx_dot.cfg" + dot_flamelet_ch4_hx.test_iter = 10 + dot_flamelet_ch4_hx.command = TestCase.Command("mpirun -n 2", "SU2_DOT_AD") + dot_flamelet_ch4_hx.timeout = 20000 + dot_flamelet_ch4_hx.reference_file = "of_grad.csv.ref" + dot_flamelet_ch4_hx.test_file = "of_grad.csv" + dot_flamelet_ch4_hx.comp_threshold = 1e-6 + dot_flamelet_ch4_hx.tol_file_percent = 0.1 + pass_list.append(dot_flamelet_ch4_hx.run_filediff()) + test_list.append(dot_flamelet_ch4_hx) + + # 2D planar laminar premixed flame on isothermal burner with conjugate heat transfer (restart) + # This test restarts on the output of test discadj_flamelet_ch4_cht and + # will only pass if test discadj_flamelet_ch4_cht passes. + dot_flamelet_ch4_cht = TestCase('dot_flamelet_ch4_cht') + dot_flamelet_ch4_cht.cfg_dir = "flamelet/04_laminar_premixed_ch4_flame_cht_ad" + dot_flamelet_ch4_cht.cfg_file = "lam_prem_ch4_cht_dot_master.cfg" + dot_flamelet_ch4_cht.multizone = True + dot_flamelet_ch4_cht.test_iter = 10 + dot_flamelet_ch4_cht.command = TestCase.Command("mpirun -n 2", "SU2_DOT_AD") + dot_flamelet_ch4_cht.timeout = 20000 + dot_flamelet_ch4_cht.reference_file = "of_grad.csv.ref" + dot_flamelet_ch4_cht.test_file = "of_grad.csv" + dot_flamelet_ch4_cht.comp_threshold = 1e-6 + dot_flamelet_ch4_cht.tol_file_percent = 0.1 + pass_list.append(dot_flamelet_ch4_cht.run_filediff()) + test_list.append(dot_flamelet_ch4_cht) + ################################################## ### Structural Adjoint - Topology Optimization ### ################################################## diff --git a/TestCases/species_transport/multizone/zone1.cfg b/TestCases/species_transport/multizone/zone1.cfg index 50c6e455e4c..c18f0565b1e 100644 --- a/TestCases/species_transport/multizone/zone1.cfg +++ b/TestCases/species_transport/multizone/zone1.cfg @@ -19,7 +19,6 @@ INLET_FILENAME= inlet_venturi.dat INC_INLET_TYPE= VELOCITY_INLET VELOCITY_INLET MARKER_INLET= ( gas_inlet1, 300.0, 1.0, 1.0, 0.0, 0.0,\ air_axial_inlet1, 305.0, 1.0, 0.0 -1.0, 0.0 ) -SPECIES_USE_STRONG_BC= NO MARKER_INLET_SPECIES= (gas_inlet1, 0.5, 0.5,\ air_axial_inlet1, 0.6, 0.0 ) MARKER_PLOTTING= (NONE) diff --git a/TestCases/species_transport/multizone/zone2.cfg b/TestCases/species_transport/multizone/zone2.cfg index 622685ab51f..2ad80c3dd58 100644 --- a/TestCases/species_transport/multizone/zone2.cfg +++ b/TestCases/species_transport/multizone/zone2.cfg @@ -18,7 +18,6 @@ SPECIFIED_INLET_PROFILE= NO INLET_FILENAME= inlet_venturi.dat INC_INLET_TYPE= VELOCITY_INLET MARKER_INLET= ( air_axial_inlet2, 310.0, 1.2, 0.0 -1.0, 0.0 ) -SPECIES_USE_STRONG_BC= NO MARKER_INLET_SPECIES= (air_axial_inlet2, 0.6, 0.0 ) INC_OUTLET_TYPE= PRESSURE_OUTLET MARKER_OUTLET= ( outlet2, 0.0 ) diff --git a/TestCases/species_transport/passive_transport_validation/passive_transport.cfg b/TestCases/species_transport/passive_transport_validation/passive_transport.cfg index 09c870157d6..a618a1e4e1c 100644 --- a/TestCases/species_transport/passive_transport_validation/passive_transport.cfg +++ b/TestCases/species_transport/passive_transport_validation/passive_transport.cfg @@ -52,7 +52,7 @@ MARKER_SYM= ( top, bottom ) INC_INLET_TYPE= VELOCITY_INLET, VELOCITY_INLET MARKER_INLET= ( gas_inlet, 300, 1, 1.0, 0.0, 0.0, air_inlet, 300, 1, 1.0, 0.0, 0.0 ) % -SPECIES_USE_STRONG_BC= YES +MARKER_SPECIES_STRONG_BC= (gas_inlet, air_inlet, outlet) MARKER_INLET_SPECIES = (gas_inlet, 1.0, air_inlet, 0.0 ) % INC_OUTLET_TYPE= PRESSURE_OUTLET diff --git a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi.cfg b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi.cfg index 6a223174a54..d39ae5957ac 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi.cfg @@ -48,7 +48,6 @@ 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 ) -SPECIES_USE_STRONG_BC= NO MARKER_INLET_SPECIES= (gas_inlet, 1.0,\ air_axial_inlet, 0.6 ) % diff --git a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_boundedscalar.cfg b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_boundedscalar.cfg index b38b3390754..4f5a9a62460 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_boundedscalar.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_boundedscalar.cfg @@ -48,7 +48,6 @@ 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 ) -SPECIES_USE_STRONG_BC= NO MARKER_INLET_SPECIES= (gas_inlet, 1.0,\ air_axial_inlet, 0.6 ) % diff --git a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel.cfg b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel.cfg index 7b022511d86..5ad7433c102 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel.cfg @@ -60,7 +60,6 @@ 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 ) -SPECIES_USE_STRONG_BC= NO MARKER_INLET_SPECIES= (gas_inlet, 1.0,\ air_axial_inlet, 0.6) % diff --git a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.cfg b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.cfg index 662646b3134..f49ea42aa8e 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_TURBULENT_MARKERS.cfg @@ -64,7 +64,7 @@ INLET_FILENAME= inlet_venturi.dat INC_INLET_TYPE= VELOCITY_INLET VELOCITY_INLET MARKER_INLET= ( gas_inlet, 300, 15.0, 1.0, 0.0, 0.0,\ air_axial_inlet, 300, 10.0, 0.0, -1.0, 0.0 ) -SPECIES_USE_STRONG_BC= YES +MARKER_SPECIES_STRONG_BC= (gas-inlet, air_axial_inlet, wall, outlet) MARKER_INLET_SPECIES= (gas_inlet, 1.0, \ air_axial_inlet, 0.0) % diff --git a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_boundedscalar.cfg b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_boundedscalar.cfg index 8a4c8d6a40d..319ab8e778e 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_boundedscalar.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_boundedscalar.cfg @@ -60,7 +60,6 @@ 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 ) -SPECIES_USE_STRONG_BC= NO MARKER_INLET_SPECIES= (gas_inlet, 1.0,\ air_axial_inlet, 0.6) % diff --git a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_heatcapacity_H2.cfg b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_heatcapacity_H2.cfg index e208b15fa59..93ad887d155 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_heatcapacity_H2.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_heatcapacity_H2.cfg @@ -70,7 +70,7 @@ INLET_FILENAME= inlet_venturi.dat INC_INLET_TYPE= VELOCITY_INLET VELOCITY_INLET MARKER_INLET= ( gas_inlet, 400, 5.0, 1.0, 0.0, 0.0,\ air_axial_inlet, 300, 5.0, 0.0, -1.0, 0.0 ) -SPECIES_USE_STRONG_BC= YES +MARKER_SPECIES_STRONG_BC= (gas_inlet, air_axial_inlet, wall, outlet) MARKER_INLET_SPECIES= (gas_inlet, 1.0, \ air_axial_inlet, 0.6) % diff --git a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.cfg b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.cfg index d941f96705f..2cd26d57e2d 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_heatcapacity_H2_ND.cfg @@ -71,7 +71,7 @@ INLET_FILENAME= inlet_venturi.dat INC_INLET_TYPE= VELOCITY_INLET VELOCITY_INLET MARKER_INLET= ( gas_inlet, 400, 5.0, 1.0, 0.0, 0.0,\ air_axial_inlet, 300, 5.0, 0.0, -1.0, 0.0 ) -SPECIES_USE_STRONG_BC= YES +MARKER_SPECIES_STRONG_BC= (gas_inlet, air_axial_inlet, wall, outlet) MARKER_INLET_SPECIES= (gas_inlet, 1.0, \ air_axial_inlet, 0.6) % diff --git a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_viscosity.cfg b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_viscosity.cfg index dd8017ec2a6..7df939a97d6 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_viscosity.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species2_primitiveVenturi_mixingmodel_viscosity.cfg @@ -65,7 +65,6 @@ INLET_FILENAME= inlet_venturi.dat INC_INLET_TYPE= VELOCITY_INLET VELOCITY_INLET MARKER_INLET= ( gas_inlet, 300, 5.0, 1.0, 0.0, 0.0,\ air_axial_inlet, 300, 3.0, 0.0, -1.0, 0.0 ) -SPECIES_USE_STRONG_BC= NO MARKER_INLET_SPECIES= (gas_inlet, 1.0, \ air_axial_inlet, 0.6) % diff --git a/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_inletFile.cfg b/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_inletFile.cfg index 6de6854f3ee..06b085e6fdf 100644 --- a/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_inletFile.cfg +++ b/TestCases/species_transport/venturi_primitive_3species/species3_primitiveVenturi_inletFile.cfg @@ -51,7 +51,6 @@ 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 ) -SPECIES_USE_STRONG_BC= NO MARKER_INLET_SPECIES= (gas_inlet, 0.5, 0.5,\ air_axial_inlet, 0.6, 0.0 ) % diff --git a/TestCases/user_defined_functions/lam_flatplate.cfg b/TestCases/user_defined_functions/lam_flatplate.cfg index c2544a6851e..457ef91fc10 100644 --- a/TestCases/user_defined_functions/lam_flatplate.cfg +++ b/TestCases/user_defined_functions/lam_flatplate.cfg @@ -37,7 +37,7 @@ CUSTOM_OUTPUTS= 'velocity : Macro{sqrt(pow(VELOCITY_X, 2) + pow(VELOCITY_Y, 2) + probe1 : Probe{$velocity}[0.005, 0.005, 0.05]' % % "COMBO" is the name and group of the output for the objective function -% (regardless of definition). "CUSTOM" is the group for all custom outpus. +% (regardless of definition). "CUSTOM" is the group for all custom outputs. SCREEN_OUTPUT= INNER_ITER, RMS_DENSITY, RMS_ENERGY, LINSOL_RESIDUAL, FORCE_Z,\ SURFACE_MASSFLOW, SURFACE_TOTAL_TEMPERATURE, avg_vel, dev_vel, probe1, COMBO HISTORY_OUTPUT = ITER, AERO_COEFF, FLOW_COEFF, FLOW_COEFF_SURF, CUSTOM, COMBO diff --git a/config_template.cfg b/config_template.cfg index a494d47eea7..e36c6d9cfc3 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -317,7 +317,7 @@ SEMI_SPAN= 0.0 % ---- NONEQUILIBRIUM GAS, IDEAL GAS, POLYTROPIC, VAN DER WAALS AND PENG ROBINSON CONSTANTS, CoolProp library -------% % % Fluid model (STANDARD_AIR, IDEAL_GAS, VW_GAS, PR_GAS, -% CONSTANT_DENSITY, INC_IDEAL_GAS, INC_IDEAL_GAS_POLY, MUTATIONPP, SU2_NONEQ, FLUID_MIXTURE, COOLPROP) +% CONSTANT_DENSITY, INC_IDEAL_GAS, INC_IDEAL_GAS_POLY, MUTATIONPP, SU2_NONEQ, FLUID_MIXTURE, COOLPROP, FLUID_FLAMELET) FLUID_MODEL= STANDARD_AIR % To find all available fluid name for CoolProp library, clikc the following link: % http://www.coolprop.org/fluid_properties/PurePseudoPure.html#list-of-fluids @@ -375,7 +375,7 @@ INLET_GAS_COMPOSITION = (0.77, 0.23, 0.0, 0.0, 0.0) % % --------------------------- VISCOSITY MODEL ---------------------------------% % -% Viscosity model (SUTHERLAND, CONSTANT_VISCOSITY, POLYNOMIAL_VISCOSITY). +% Viscosity model (SUTHERLAND, CONSTANT_VISCOSITY, POLYNOMIAL_VISCOSITY, FLAMELET). VISCOSITY_MODEL= SUTHERLAND % % Molecular Viscosity that would be constant (1.716E-5 by default) @@ -397,7 +397,7 @@ MU_POLYCOEFFS= (0.0, 0.0, 0.0, 0.0, 0.0) % --------------------------- THERMAL CONDUCTIVITY MODEL ----------------------% % % Laminar Conductivity model (CONSTANT_CONDUCTIVITY, CONSTANT_PRANDTL, -% POLYNOMIAL_CONDUCTIVITY). +% POLYNOMIAL_CONDUCTIVITY, FLAMELET). CONDUCTIVITY_MODEL= CONSTANT_PRANDTL % % Molecular Thermal Conductivity that would be constant (0.0257 by default) @@ -707,10 +707,10 @@ TIME_DISCRE_RADIATION = EULER_IMPLICIT % --------------------- SPECIES TRANSPORT SIMULATION --------------------------% % -% Specify scalar transport model (NONE, SPECIES_TRANSPORT) +% Specify scalar transport model (NONE, SPECIES_TRANSPORT, FLAMELET) KIND_SCALAR_MODEL= NONE % -% Mass diffusivity model (CONSTANT_DIFFUSIVITY) +% Mass diffusivity model (CONSTANT_DIFFUSIVITY, FLAMELET) DIFFUSIVITY_MODEL= CONSTANT_DIFFUSIVITY % % Mass diffusivity if DIFFUSIVITY_MODEL= CONSTANT_DIFFUSIVITY is chosen. D_air ~= 0.001 @@ -724,9 +724,6 @@ SCHMIDT_NUMBER_TURBULENT= 0.7 % For N species, N-1 transport equations are solved, the last one Y_N is solved algebraically as 1-(sum of the species 1 to (N-1)) MARKER_INLET_SPECIES= (inlet, 0.5, ..., inlet2, 0.6, ...) % -% Use strong inlet and outlet BC in the species solver -SPECIES_USE_STRONG_BC= NO -% % Convective numerical method for species transport (SCALAR_UPWIND, BOUNDED_SCALAR) CONV_NUM_METHOD_SPECIES= SCALAR_UPWIND % @@ -755,6 +752,35 @@ SPECIES_CLIPPING_MAX= 1.0, ... % Minimum values for scalar clipping SPECIES_CLIPPING_MIN= 0.0, ... +% --------------------- FLAMELET MODEL -----------------------------% +% +% Names of the user defined (auxiliary) transport equations. +USER_SCALAR_NAMES= (Y-CO) +% +% Names of the source terms for the user defined transport equations. The source term +% has the form +% S = S_production + S_consumption * Y, with Y the mass fraction of the user scalar +% S = S_production + S_consumption * X, with X the mole fraction of the user scalar +% The second term can have the name "zero" or "null" to indicate that only a total source +% will be used. +USER_SOURCE_NAMES= ( ProdRateTot_Y-CO, zero) +% +% Names of variables from the LUT that will be looked up for visualization +LOOKUP_NAMES= (MolarWeightMix, Conductivity, HeatRelease, Diffusivity) +% +% filename of the lookup table +FILENAME_LUT= fgm_ch4.drg +% +% flamelet initialization +% the flame is initialized using a plane, defined by a point and a normal. On one side, the solution is initialized +% using 'burnt' conditions and on the other side 'unburnt' conditions. The normal points in the direction of the 'burnt' +% condition. +% (x1,x2,x3) = point on the plane. +% (x4,x5,x6) = normal of the plane. +% (x7) = thickness of the reaction zone, this is the transition from unburnt to burnt conditions. +% (x8) = Thickness of the 'burnt' zone, after this length, the conditions will be 'unburnt' again. +FLAME_INIT= (0.004, 0.0, 0.0, 1.0, 0.0, 0.0, 0.2e-3, 1.0) +% % --------------------- INVERSE DESIGN SIMULATION -----------------------------% % % Evaluate an inverse design problem using Cp (NO, YES) @@ -991,6 +1017,9 @@ CATALYTIC_WALL= ( NONE ) % SA Model: (inlet_marker1, NuFactor1, inlet_marker2, NuFactor2, ...) % SST Model: (inlet_marker1, TurbIntensity1, RatioTurbLamViscosity1, inlet_marker2, TurbIntensity2, RatioTurbLamViscosity2, ...) MARKER_INLET_TURBULENT= (inlet1, 0.05, 15, inlet2, 0.02, ...) +% +% list of markers species transport and flamelet model where strong boundary conditions should be used +MARKER_SPECIES_STRONG_BC= (inlet, wall) % ------------------------ WALL ROUGHNESS DEFINITION --------------------------% % The equivalent sand grain roughness height (k_s) on each of the wall. This must be in m.