Merge branch 'develop' into preacc_performance

Author: jblueh

Common/include/CConfig.hpp

@@ -909,6 +909,7 @@ class CConfig {
   Tke_FreeStreamND,           /*!< \brief Farfield kinetic energy (external flow). */
   Omega_FreeStreamND,         /*!< \brief Specific dissipation (external flow). */
   Omega_FreeStream;           /*!< \brief Specific dissipation (external flow). */
+  bool Variable_Density;      /*!< \brief Variable density for incompressible flow. */
   unsigned short nElectric_Constant;    /*!< \brief Number of different electric constants. */
   su2double *Electric_Constant;         /*!< \brief Dielectric constant modulus. */
   su2double Knowles_B,                  /*!< \brief Knowles material model constant B. */
@@ -1229,9 +1230,10 @@ class CConfig {
   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. */
 
+  string* controlling_variable_names;
+  string* cv_source_names;
   vector<string> 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* lookup_names;         /*!< \brief Names of passive look-up variables. */
   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. */
 
@@ -2151,6 +2153,19 @@ class CConfig {
    */
   unsigned short GetNUserScalars() const { return n_user_scalars; }
 
+  /*!
+   * \brief Get the name of a specific controlling variable.
+   */
+  const string& GetControllingVariableName(unsigned short i_cv) const {
+    return controlling_variable_names[i_cv];
+  }
+
+  /*!
+   * \brief Get the name of the source term variable for a specific controlling variable.
+   */
+  const string& GetControllingVariableSourceName(unsigned short i_cv) const {
+    return cv_source_names[i_cv];
+  }
   /*!
    * \brief Get the name of the user scalar.
    */
@@ -2175,13 +2190,7 @@ class CConfig {
   /*!
    * \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; }
+  const string& GetLookupName(unsigned short i_lookup) const { return lookup_names[i_lookup]; }
 
   /*!
    * \brief Get the Young's modulus of elasticity.
@@ -3884,6 +3893,12 @@ class CConfig {
    */
   INC_DENSITYMODEL GetKind_DensityModel() const { return Kind_DensityModel; }
 
+  /*!
+   * \brief Selection of variable density option for incompressible flows.
+   * \return Flag for variable density for incompressible flows.
+   */
+  bool GetVariable_Density_Model() const { return Variable_Density; }
+
   /*!
    * \brief Flag for whether to solve the energy equation for incompressible flows.
    * \return Flag for energy equation

Common/include/containers/CLookUpTable.hpp

@@ -76,6 +76,8 @@ class CLookUpTable {
    */
   su2vector<su2activematrix> table_data;
 
+  double memory_footprint_data = 0; /*!< \brief Memory footprint of the loaded table data. */
+
   /*! \brief
    * Holds all connectivity data stored in the table for each level. First index
    * addresses the variable while second index addresses the point.
@@ -216,7 +218,7 @@ class CLookUpTable {
    */
   void Linear_Interpolation(const su2double val_CV3, const unsigned long lower_level, const unsigned long upper_level,
                             std::vector<su2double>& lower_values, std::vector<su2double>& upper_values,
-                            std::vector<su2double*>& var_vals) const;
+                            std::vector<su2double>& 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).
@@ -345,7 +347,7 @@ class CLookUpTable {
    * \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<std::string>& val_names_var, std::vector<su2double*>& val_vars,
+  unsigned long LookUp_XYZ(const std::vector<std::string>& val_names_var, std::vector<su2double>& val_vars,
                            su2double val_CV1, su2double val_CV2, su2double val_CV3 = 0);
 
   /*!

Common/include/containers/CTrapezoidalMap.hpp

@@ -54,12 +54,14 @@ class CTrapezoidalMap {
    * same band. Used to sort the edges */
   su2vector<std::vector<std::pair<su2double, unsigned long> > > y_edge_at_band_mid;
 
+  double memory_footprint = 0;
+
  public:
   CTrapezoidalMap() = default;
 
   CTrapezoidalMap(const su2double* samples_x, const su2double* samples_y, const unsigned long size,
                   const std::vector<std::array<unsigned long, 2> >& edges,
-                  const su2vector<std::vector<unsigned long> >& edge_to_triangle);
+                  const su2vector<std::vector<unsigned long> >& edge_to_triangle, bool display = false);
 
   /*!
    * \brief return the index to the triangle that contains the coordinates (val_x,val_y)
@@ -97,4 +99,10 @@ class CTrapezoidalMap {
   inline bool IsInsideHullX(su2double val_x) {
     return (val_x >= unique_bands_x.front()) && (val_x <= unique_bands_x.back());
   }
+
+  /*!
+   * \brief get memory footprint of trapezoidal map.
+   * \return - memory footprint in mega bytes.
+   */
+  double GetMemoryFootprint() const { return memory_footprint; }
 };

Common/include/linear_algebra/CSysSolve.hpp

@@ -168,6 +168,7 @@ class CSysSolve {
    * \brief Modified Gram-Schmidt orthogonalization
    * \author Based on Kesheng John Wu's mgsro subroutine in Saad's SPARSKIT
    *
+   * \param[in] shared_hsbg - if the Hessenberg matrix is shared by multiple threads
    * \param[in] i - index indicating which vector in w is being orthogonalized
    * \param[in,out] Hsbg - the upper Hessenberg begin updated
    * \param[in,out] w - the (i+1)th vector of w is orthogonalized against the
@@ -181,7 +182,7 @@ class CSysSolve {
    * vector is kept in nrm0 and updated after operating with each vector
    *
    */
-  void ModGramSchmidt(int i, su2matrix<ScalarType>& Hsbg, std::vector<VectorType>& w) const;
+  void ModGramSchmidt(bool shared_hsbg, int i, su2matrix<ScalarType>& Hsbg, std::vector<VectorType>& w) const;
 
   /*!
    * \brief writes header information for a CSysSolve residual history

Common/include/linear_algebra/vector_expressions.hpp

@@ -146,28 +146,30 @@ MAKE_UNARY_FUN(sign, sign_, sign_impl)
 
 /*--- Macro to create expressions and overloads for binary functions. ---*/
 
-#define MAKE_BINARY_FUN(FUN, EXPR, IMPL)                                                                             \
-  /*!--- Expression class. ---*/                                                                                     \
-  template <class U, class V, class Scalar>                                                                          \
-  class EXPR : public CVecExpr<EXPR<U, V, Scalar>, Scalar> {                                                         \
-    store_t<const U> u;                                                                                              \
-    store_t<const V> v;                                                                                              \
-                                                                                                                     \
-   public:                                                                                                           \
-    static constexpr bool StoreAsRef = false;                                                                        \
-    FORCEINLINE EXPR(const U& u_, const V& v_) : u(u_), v(v_) {}                                                     \
-    FORCEINLINE auto operator[](size_t i) const RETURNS(IMPL(u[i], v[i]))                                            \
-  };                                                                                                                 \
-  /*!--- Vector with vector function overload. ---*/                                                                 \
-  template <class U, class V, class S>                                                                               \
-  FORCEINLINE auto FUN(const CVecExpr<U, S>& u, const CVecExpr<V, S>& v)                                             \
-      RETURNS(EXPR<U, V, S>(u.derived(), v.derived())) /*!--- Vector with scalar function overload. ---*/            \
-      template <class U, class S>                                                                                    \
-      FORCEINLINE auto FUN(const CVecExpr<U, S>& u, decay_t<S> v)                                                    \
-          RETURNS(EXPR<U, Bcast<S>, S>(u.derived(), Bcast<S>(v))) /*!--- Scalar with vector function overload. ---*/ \
-      template <class S, class V>                                                                                    \
-      FORCEINLINE auto FUN(decay_t<S> u, const CVecExpr<V, S>& v)                                                    \
-          RETURNS(EXPR<Bcast<S>, V, S>(Bcast<S>(u), v.derived()))
+// clang-format off
+#define MAKE_BINARY_FUN(FUN, EXPR, IMPL)                                                                              \
+  /*!--- Expression class. ---*/                                                                                      \
+  template <class U, class V, class Scalar>                                                                           \
+  class EXPR : public CVecExpr<EXPR<U, V, Scalar>, Scalar> {                                                          \
+    store_t<const U> u;                                                                                               \
+    store_t<const V> v;                                                                                               \
+                                                                                                                      \
+   public:                                                                                                            \
+    static constexpr bool StoreAsRef = false;                                                                         \
+    FORCEINLINE EXPR(const U& u_, const V& v_) : u(u_), v(v_) {}                                                      \
+    FORCEINLINE auto operator[](size_t i) const RETURNS(IMPL(u[i], v[i]))                                             \
+  };                                                                                                                  \
+  /*!--- Vector with vector function overload. ---*/                                                                  \
+  template <class U, class V, class S>                                                                                \
+  FORCEINLINE auto FUN(const CVecExpr<U, S>& u, const CVecExpr<V, S>& v)                                              \
+      RETURNS(EXPR<U, V, S>(u.derived(), v.derived()))                                                                \
+  /*!--- Vector with scalar function overload. ---*/                                                                  \
+  template <class U, class S>                                                                                         \
+  FORCEINLINE auto FUN(const CVecExpr<U, S>& u, decay_t<S> v) RETURNS(EXPR<U, Bcast<S>, S>(u.derived(), Bcast<S>(v))) \
+  /*!--- Scalar with vector function overload. ---*/                                                                  \
+  template <class S, class V>                                                                                         \
+  FORCEINLINE auto FUN(decay_t<S> u, const CVecExpr<V, S>& v) RETURNS(EXPR<Bcast<S>, V, S>(Bcast<S>(u), v.derived()))
+// clang-format on
 
 /*--- std::max/min have issues (because they return by reference).
  * fmin and fmax return by value and thus are fine, but they would force

Common/include/option_structure.hpp

@@ -628,11 +628,13 @@ enum class INC_DENSITYMODEL {
   CONSTANT,   /*!< \brief Constant density. */
   BOUSSINESQ, /*!< \brief Boussinesq density model. */
   VARIABLE,   /*!< \brief Variable density model. */
+  FLAMELET,   /*!< \brief Density according to flamelet manifold. */
 };
 static const MapType<std::string, INC_DENSITYMODEL> DensityModel_Map = {
   MakePair("CONSTANT", INC_DENSITYMODEL::CONSTANT)
   MakePair("BOUSSINESQ", INC_DENSITYMODEL::BOUSSINESQ)
   MakePair("VARIABLE", INC_DENSITYMODEL::VARIABLE)
+  MakePair("FLAMELET", INC_DENSITYMODEL::FLAMELET)
 };
 
 /*!
@@ -1317,6 +1319,7 @@ static const MapType<std::string, SPECIES_MODEL> Species_Model_Map = {
 enum FLAMELET_SCALAR_VARIABLES {
   I_PROGVAR,
   I_ENTH,
+  I_MIXFRAC,
 };
 
 /*!
@@ -1326,6 +1329,15 @@ enum FLAMELET_SCALAR_SOURCES {
   I_SRC_TOT_PROGVAR
 };
 
+/*!
+ * \brief Look-up operations for the flamelet scalar solver.
+ */
+enum FLAMELET_LOOKUP_OPS {
+  TD,       /*!< \brief Thermochemical properties (temperature, density, diffusivity, etc.). */
+  SOURCES,  /*!< \brief Scalar source terms (controlling variables, passive species).*/
+  LOOKUP,   /*!< \brief Passive look-up variables specified in config. */
+};
+
 /*!
  * \brief Types of subgrid scale models
  */

Common/src/CConfig.cpp

@@ -2101,11 +2101,14 @@ 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"));
+  /* !\brief CONTROLLING_VARIABLE_NAMES \n DESCRIPTION: Names of the variables used as inputs for the data regression method in flamelet or data-driven fluid models. */
+  addStringListOption("CONTROLLING_VARIABLE_NAMES", n_control_vars, controlling_variable_names);
+
+  /* !\brief CONTROLLING_VARIABLE_SOURCE_NAMES \n DESCRIPTION: Names of the variables in the flamelet manifold corresponding to the source terms of the controlling variables. */
+  addStringListOption("CONTROLLING_VARIABLE_SOURCE_NAMES", n_control_vars, cv_source_names);
 
   /* DESCRIPTION: Names of the passive lookup variables for flamelet LUT */
-  addStringListOption("LOOKUP_NAMES", n_lookups, table_lookup_names);
+  addStringListOption("LOOKUP_NAMES", n_lookups, lookup_names);
 
   /* DESCRIPTION: Names of the user transport equations solved in the flamelet problem. */
   addStringListOption("USER_SCALAR_NAMES", n_user_scalars, user_scalar_names);
@@ -3840,6 +3843,8 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i
   SetDefaultIfEmpty(Prandtl_Turb, nPrandtl_Turb, Prandtl_Turb_Default);
   SetDefaultIfEmpty(Constant_Lewis_Number, nConstant_Lewis_Number, Lewis_Number_Default);
 
+  Variable_Density = ((Kind_DensityModel == INC_DENSITYMODEL::VARIABLE) || (Kind_DensityModel == INC_DENSITYMODEL::FLAMELET));
+
   /*--- Check whether inputs for FLUID_MIXTURE are correctly specified. ---*/
 
     if (Kind_FluidModel == FLUID_MIXTURE) {
@@ -3935,8 +3940,8 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i
                        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",
+      if (!Variable_Density) {
+        SU2_MPI::Error("The use of FLUID_FLAMELET requires the INC_DENSITY_MODEL option to be VARIABLE or FLAMELET",
                        CURRENT_FUNCTION);
       }
 
@@ -5480,8 +5485,10 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i
 
   /*--- 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;
+    /*--- The controlling variables are progress variable, total enthalpy, and optionally mixture fraction ---*/
+    //n_control_vars = nSpecies - n_user_scalars;
+    if (n_control_vars != (nSpecies - n_user_scalars))
+      SU2_MPI::Error("Number of initial species incompatbile with number of controlling variables and user scalars.", CURRENT_FUNCTION);
     /*--- We can have additional user defined transported scalars ---*/
     n_scalars = n_control_vars + n_user_scalars;
   }