Merge pull request #1010 from su2code/thermophysical

Divides and Moves the Fluid Property Classes

Author: economon

Common/include/CConfig.hpp

@@ -783,8 +783,6 @@ class CConfig {
   unsigned short
   Console_Output_Verb,  /*!< \brief Level of verbosity for console output */
   Kind_Average;         /*!< \brief Particular average for the marker analyze. */
-  unsigned short
-  nPolyCoeffs;          /*!< \brief Number of coefficients in temperature polynomial fits for fluid models. */
   su2double Gamma,      /*!< \brief Ratio of specific heats of the gas. */
   Bulk_Modulus,         /*!< \brief Value of the bulk modulus for incompressible flows. */
   Beta_Factor,          /*!< \brief Value of the epsilon^2 multiplier for Beta for the incompressible preconditioner. */
@@ -818,14 +816,14 @@ class CConfig {
   Mu_Temperature_Ref,    /*!< \brief Reference temperature for Sutherland model.  */
   Mu_Temperature_RefND,  /*!< \brief Non-dimensional reference temperature for Sutherland model.  */
   Mu_S,                  /*!< \brief Reference S for Sutherland model.  */
-  Mu_SND,                /*!< \brief Non-dimensional reference S for Sutherland model.  */
-  *CpPolyCoefficients,     /*!< \brief Definition of the temperature polynomial coefficients for specific heat Cp. */
-  *MuPolyCoefficients,     /*!< \brief Definition of the temperature polynomial coefficients for viscosity. */
-  *KtPolyCoefficients,     /*!< \brief Definition of the temperature polynomial coefficients for thermal conductivity. */
-  *CpPolyCoefficientsND,   /*!< \brief Definition of the non-dimensional temperature polynomial coefficients for specific heat Cp. */
-  *MuPolyCoefficientsND,   /*!< \brief Definition of the non-dimensional temperature polynomial coefficients for viscosity. */
-  *KtPolyCoefficientsND,   /*!< \brief Definition of the non-dimensional temperature polynomial coefficients for thermal conductivity. */
-  Thermal_Conductivity_Solid,      /*!< \brief Thermal conductivity in solids. */
+  Mu_SND;                /*!< \brief Non-dimensional reference S for Sutherland model.  */
+  su2double* CpPolyCoefficients;     /*!< \brief Definition of the temperature polynomial coefficients for specific heat Cp. */
+  su2double* MuPolyCoefficients;     /*!< \brief Definition of the temperature polynomial coefficients for viscosity. */
+  su2double* KtPolyCoefficients;     /*!< \brief Definition of the temperature polynomial coefficients for thermal conductivity. */
+  array<su2double, N_POLY_COEFFS> CpPolyCoefficientsND{0.0};   /*!< \brief Definition of the non-dimensional temperature polynomial coefficients for specific heat Cp. */
+  array<su2double, N_POLY_COEFFS>MuPolyCoefficientsND{0.0};   /*!< \brief Definition of the non-dimensional temperature polynomial coefficients for viscosity. */
+  array<su2double, N_POLY_COEFFS>KtPolyCoefficientsND{0.0};   /*!< \brief Definition of the non-dimensional temperature polynomial coefficients for thermal conductivity. */
+  su2double Thermal_Conductivity_Solid,      /*!< \brief Thermal conductivity in solids. */
   Thermal_Diffusivity_Solid,       /*!< \brief Thermal diffusivity in solids. */
   Temperature_Freestream_Solid,    /*!< \brief Temperature in solids at freestream conditions. */
   Density_Solid,                   /*!< \brief Total density in solids. */
@@ -1018,9 +1016,9 @@ class CConfig {
   su2double FinalRotation_Rate_Z;       /*!< \brief Final rotation rate Z if Ramp rotating frame is activated. */
   su2double FinalOutletPressure;        /*!< \brief Final outlet pressure if Ramp outlet pressure is activated. */
   su2double MonitorOutletPressure;      /*!< \brief Monitor outlet pressure if Ramp outlet pressure is activated. */
-  su2double *default_cp_polycoeffs;     /*!< \brief Array for specific heat polynomial coefficients. */
-  su2double *default_mu_polycoeffs;     /*!< \brief Array for viscosity polynomial coefficients. */
-  su2double *default_kt_polycoeffs;     /*!< \brief Array for thermal conductivity polynomial coefficients. */
+  array<su2double, N_POLY_COEFFS> default_cp_polycoeffs{0.0};     /*!< \brief Array for specific heat polynomial coefficients. */
+  array<su2double, N_POLY_COEFFS> default_mu_polycoeffs{0.0};     /*!< \brief Array for viscosity polynomial coefficients. */
+  array<su2double, N_POLY_COEFFS> default_kt_polycoeffs{0.0};     /*!< \brief Array for thermal conductivity polynomial coefficients. */
   su2double *ExtraRelFacGiles;          /*!< \brief coefficient for extra relaxation factor for Giles BC*/
   bool Body_Force;                      /*!< \brief Flag to know if a body force is included in the formulation. */
   su2double *Body_Force_Vector;         /*!< \brief Values of the prescribed body force vector. */
@@ -3818,7 +3816,7 @@ class CConfig {
    * \brief Get the number of coefficients in the temperature polynomial models.
    * \return The the number of coefficients in the temperature polynomial models.
    */
-  unsigned short GetnPolyCoeffs(void) const { return nPolyCoeffs; }
+  unsigned short GetnPolyCoeffs(void) const { return N_POLY_COEFFS; }
 
   /*!
    * \brief Get the temperature polynomial coefficient for specific heat Cp.
@@ -3852,7 +3850,7 @@ class CConfig {
    * \brief Get the temperature polynomial coefficients for viscosity.
    * \return Non-dimensional temperature polynomial coefficients for viscosity.
    */
-  su2double* GetMu_PolyCoeffND(void) { return MuPolyCoefficientsND; }
+  const su2double* GetMu_PolyCoeffND(void) const { return MuPolyCoefficientsND.data(); }
 
   /*!
    * \brief Get the temperature polynomial coefficient for thermal conductivity.
@@ -3872,7 +3870,7 @@ class CConfig {
    * \brief Get the temperature polynomial coefficients for thermal conductivity.
    * \return Non-dimensional temperature polynomial coefficients for thermal conductivity.
    */
-  su2double* GetKt_PolyCoeffND(void) { return KtPolyCoefficientsND; }
+  const su2double* GetKt_PolyCoeffND(void) const { return KtPolyCoefficientsND.data(); }
 
   /*!
    * \brief Set the value of the non-dimensional constant viscosity.

Common/include/fem/fem_geometry_structure.hpp

@@ -174,7 +174,7 @@ class CSortFaces {
  * \author E. van der Weide
  * \version 7.0.5 "Blackbird"
  */
-class CSurfaceElementFEM;   // Forward declaration to avoid problems.
+struct CSurfaceElementFEM;   // Forward declaration to avoid problems.
 struct CSortBoundaryFaces {
  /*!
   * \brief Operator used for the comparison.

Common/include/option_structure.hpp

@@ -136,6 +136,12 @@ const int SU2_CONN_SKIP   = 2;   /*!< \brief Offset to skip the globalID and VTK
 
 const su2double COLORING_EFF_THRESH = 0.875;  /*!< \brief Below this value fallback strategies are used instead. */
 
+/*--- All temperature polynomial fits for the fluid models currently
+   assume a quartic form (5 coefficients). For example,
+   Cp(T) = b0 + b1*T + b2*T^2 + b3*T^3 + b4*T^4. By default, all coeffs
+   are set to zero and will be properly non-dim. in the solver. ---*/
+constexpr int N_POLY_COEFFS = 5; /*!< \brief Number of coefficients in temperature polynomial fits for fluid models. */;
+
 /*!
  * \brief Boolean answers
  */

Common/src/CConfig.cpp

@@ -978,13 +978,6 @@ void CConfig::SetPointersNull(void) {
 
   /*--- Initialize some default arrays to NULL. ---*/
 
-  default_cp_polycoeffs = nullptr;
-  default_mu_polycoeffs = nullptr;
-  default_kt_polycoeffs = nullptr;
-  CpPolyCoefficientsND  = nullptr;
-  MuPolyCoefficientsND  = nullptr;
-  KtPolyCoefficientsND  = nullptr;
-
   Riemann_FlowDir       = nullptr;
   Giles_FlowDir         = nullptr;
   CoordFFDBox           = nullptr;
@@ -1076,23 +1069,6 @@ void CConfig::SetRunTime_Options(void) {
 
 void CConfig::SetConfig_Options() {
 
-
-  /*--- Allocate some default arrays needed for lists of doubles. ---*/
-
-
-  /*--- All temperature polynomial fits for the fluid models currently
-   assume a quartic form (5 coefficients). For example,
-   Cp(T) = b0 + b1*T + b2*T^2 + b3*T^3 + b4*T^4. By default, all coeffs
-   are set to zero and will be properly non-dim. in the solver. ---*/
-
-  nPolyCoeffs = 5;
-  default_cp_polycoeffs = new su2double[nPolyCoeffs]();
-  default_mu_polycoeffs = new su2double[nPolyCoeffs]();
-  default_kt_polycoeffs = new su2double[nPolyCoeffs]();
-  CpPolyCoefficientsND  = new su2double[nPolyCoeffs]();
-  MuPolyCoefficientsND  = new su2double[nPolyCoeffs]();
-  KtPolyCoefficientsND  = new su2double[nPolyCoeffs]();
-
   // This config file is parsed by a number of programs to make it easy to write SU2
   // wrapper scripts (in python, go, etc.) so please do
   // the best you can to follow the established format. It's very hard to parse c++ code
@@ -1219,11 +1195,11 @@ void CConfig::SetConfig_Options() {
   /*--- Options related to temperature polynomial coefficients for fluid models. ---*/
 
   /* DESCRIPTION: Definition of the temperature polynomial coefficients for specific heat Cp. */
-  addDoubleArrayOption("CP_POLYCOEFFS", nPolyCoeffs, CpPolyCoefficients, default_cp_polycoeffs);
+  addDoubleArrayOption("CP_POLYCOEFFS", N_POLY_COEFFS, CpPolyCoefficients, default_cp_polycoeffs.data());
   /* DESCRIPTION: Definition of the temperature polynomial coefficients for specific heat Cp. */
-  addDoubleArrayOption("MU_POLYCOEFFS", nPolyCoeffs, MuPolyCoefficients, default_mu_polycoeffs);
+  addDoubleArrayOption("MU_POLYCOEFFS", N_POLY_COEFFS, MuPolyCoefficients, default_mu_polycoeffs.data());
   /* DESCRIPTION: Definition of the temperature polynomial coefficients for specific heat Cp. */
-  addDoubleArrayOption("KT_POLYCOEFFS", nPolyCoeffs, KtPolyCoefficients, default_kt_polycoeffs);
+  addDoubleArrayOption("KT_POLYCOEFFS", N_POLY_COEFFS, KtPolyCoefficients, default_kt_polycoeffs.data());
 
   /*!\brief REYNOLDS_NUMBER \n DESCRIPTION: Reynolds number (non-dimensional, based on the free-stream values). Needed for viscous solvers. For incompressible solvers the Reynolds length will always be 1.0 \n DEFAULT: 0.0 \ingroup Config */
   addDoubleOption("REYNOLDS_NUMBER", Reynolds, 0.0);
@@ -1609,9 +1585,9 @@ void CConfig::SetConfig_Options() {
   addBoolOption("CFL_ADAPT", CFL_Adapt, false);
   /* !\brief CFL_ADAPT_PARAM
    * DESCRIPTION: Parameters of the adaptive CFL number (factor down, factor up, CFL limit (min and max) )
-   * Factor down generally >1.0, factor up generally < 1.0 to cause the CFL to increase when residual is decreasing,
-   * and decrease when the residual is increasing or stalled. \ingroup Config*/
-  default_cfl_adapt[0] = 0.0; default_cfl_adapt[1] = 0.0; default_cfl_adapt[2] = 1.0; default_cfl_adapt[3] = 100.0;
+   * Factor down generally <1.0, factor up generally > 1.0 to cause the CFL to increase when the under-relaxation parameter is 1.0 
+   * and to decrease when the under-relaxation parameter is less than 0.1. Factor is multiplicative. \ingroup Config*/
+  default_cfl_adapt[0] = 1.0; default_cfl_adapt[1] = 1.0; default_cfl_adapt[2] = 10.0; default_cfl_adapt[3] = 100.0;
   addDoubleArrayOption("CFL_ADAPT_PARAM", 4, CFL_AdaptParam, default_cfl_adapt);
   /* DESCRIPTION: Reduction factor of the CFL coefficient in the adjoint problem */
   addDoubleOption("CFL_REDUCTION_ADJFLOW", CFLRedCoeff_AdjFlow, 0.8);
@@ -4549,28 +4525,28 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_
 
   if ((Kind_Solver == INC_EULER || Kind_Solver == INC_NAVIER_STOKES || Kind_Solver == INC_RANS) && (Kind_FluidModel == INC_IDEAL_GAS_POLY)) {
     su2double sum = 0.0;
-    for (unsigned short iVar = 0; iVar < nPolyCoeffs; iVar++) {
+    for (unsigned short iVar = 0; iVar < N_POLY_COEFFS; iVar++) {
       sum += GetCp_PolyCoeff(iVar);
     }
-    if ((nPolyCoeffs < 1) || (sum == 0.0))
+    if ((N_POLY_COEFFS < 1) || (sum == 0.0))
       SU2_MPI::Error(string("CP_POLYCOEFFS not set for fluid model INC_IDEAL_GAS_POLY. \n"), CURRENT_FUNCTION);
   }
 
   if (((Kind_Solver == INC_EULER || Kind_Solver == INC_NAVIER_STOKES || Kind_Solver == INC_RANS)) && (Kind_ViscosityModel == POLYNOMIAL_VISCOSITY)) {
     su2double sum = 0.0;
-    for (unsigned short iVar = 0; iVar < nPolyCoeffs; iVar++) {
+    for (unsigned short iVar = 0; iVar < N_POLY_COEFFS; iVar++) {
       sum += GetMu_PolyCoeff(iVar);
     }
-    if ((nPolyCoeffs < 1) || (sum == 0.0))
+    if ((N_POLY_COEFFS < 1) || (sum == 0.0))
       SU2_MPI::Error(string("MU_POLYCOEFFS not set for viscosity model POLYNOMIAL_VISCOSITY. \n"), CURRENT_FUNCTION);
   }
 
   if ((Kind_Solver == INC_EULER || Kind_Solver == INC_NAVIER_STOKES || Kind_Solver == INC_RANS) && (Kind_ConductivityModel == POLYNOMIAL_CONDUCTIVITY)) {
     su2double sum = 0.0;
-    for (unsigned short iVar = 0; iVar < nPolyCoeffs; iVar++) {
+    for (unsigned short iVar = 0; iVar < N_POLY_COEFFS; iVar++) {
       sum += GetKt_PolyCoeff(iVar);
     }
-    if ((nPolyCoeffs < 1) || (sum == 0.0))
+    if ((N_POLY_COEFFS < 1) || (sum == 0.0))
       SU2_MPI::Error(string("KT_POLYCOEFFS not set for conductivity model POLYNOMIAL_CONDUCTIVITY. \n"), CURRENT_FUNCTION);
   }
 
@@ -4894,6 +4870,20 @@ void CConfig::SetPostprocessing(unsigned short val_software, unsigned short val_
                    CURRENT_FUNCTION);
   }
 
+  /* Protect against using incorrect CFL adaption parameters. */
+  
+  if (CFL_Adapt && (CFL_AdaptParam[0] > 1.0)) {
+    SU2_MPI::Error(string("CFL adaption factor down should be less than 1.0."), CURRENT_FUNCTION);
+  }
+  
+  if (CFL_Adapt && (CFL_AdaptParam[1] < 1.0)) {
+    SU2_MPI::Error(string("CFL adaption factor up should be greater than 1.0."), CURRENT_FUNCTION);
+  }
+
+  if (CFL_Adapt && (CFL_AdaptParam[2] > CFL_AdaptParam[3])) {
+    SU2_MPI::Error(string("CFL adaption minimum CFL is larger than the maximum CFL."), CURRENT_FUNCTION);
+  }
+
   /*--- 0 in the config file means "disable" which can be done using a very large group. ---*/
   if (edgeColorGroupSize==0) edgeColorGroupSize = 1<<30;
 
@@ -7629,13 +7619,6 @@ CConfig::~CConfig(void) {
 
   /*--- Delete some arrays needed just for initializing options. ---*/
 
-  delete [] default_cp_polycoeffs;
-  delete [] default_mu_polycoeffs;
-  delete [] default_kt_polycoeffs;
-  delete [] CpPolyCoefficientsND;
-  delete [] MuPolyCoefficientsND;
-  delete [] KtPolyCoefficientsND;
-
   delete [] FFDTag;
   delete [] nDV_Value;
   delete [] TagFFDBox;

Common/src/wall_model.cpp

@@ -27,7 +27,7 @@
  */
 
 #include "../include/wall_model.hpp"
-#include "../../SU2_CFD/include/fluid_model.hpp"
+#include "../../SU2_CFD/include/fluid/CFluidModel.hpp"
 
 /* Prototypes for Lapack functions, if MKL or LAPACK is used. */
 #if defined (HAVE_MKL) || defined(HAVE_LAPACK)

SU2_CFD/include/fluid/CConductivityModel.hpp

@@ -0,0 +1,71 @@
+/*!
+ * \file CConductivityModel.hpp
+ * \brief Defines an interface class for thermal conductivity models.
+ * \author S. Vitale, M. Pini, G. Gori, A. Guardone, P. Colonna, T. Economon
+ * \version 7.0.5 "Blackbird"
+ *
+ * SU2 Project Website: https://su2code.github.io
+ *
+ * The SU2 Project is maintained by the SU2 Foundation
+ * (http://su2foundation.org)
+ *
+ * Copyright 2012-2020, 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 <http://www.gnu.org/licenses/>.
+ */
+
+#pragma once
+
+#include "../../../Common/include/basic_types/datatype_structure.hpp"
+
+using namespace std;
+
+/*!
+ * \class CConductivityModel
+ * \brief Interface class for defining the thermal conductivity model.
+ * \author S. Vitale, M. Pini
+ */
+class CConductivityModel {
+ public:
+  CConductivityModel() = default;
+  CConductivityModel(const CConductivityModel&) = delete;
+  void operator=(const CConductivityModel&) = delete;
+  virtual ~CConductivityModel() {}
+
+  /*!
+   * \brief return conductivity value.
+   */
+  virtual su2double GetConductivity(void) const = 0;
+
+  /*!
+   * \brief return conductivity partial derivative value.
+   */
+  virtual su2double Getdktdrho_T(void) const = 0;
+
+  /*!
+   * \brief return viscosity partial derivative value.
+   */
+  virtual su2double GetdktdT_rho(void) const = 0;
+
+  /*!
+   * \brief Set thermal conductivity.
+   */
+  virtual void SetConductivity(su2double t, su2double rho, su2double mu_lam, su2double mu_turb, su2double cp) = 0;
+
+  /*!
+   * \brief Set thermal conductivity derivatives.
+   */
+  virtual void SetDerConductivity(su2double t, su2double rho, su2double dmudrho_t, su2double dmudt_rho,
+                                  su2double cp) = 0;
+};

SU2_CFD/include/fluid/CConstantConductivity.hpp

@@ -0,0 +1,74 @@
+/*!
+ * \file CConstantConductivity.hpp
+ * \brief Defines a constant laminar thermal conductivity model.
+ * \author S. Vitale, M. Pini, G. Gori, A. Guardone, P. Colonna, T. Economon
+ * \version 7.0.5 "Blackbird"
+ *
+ * SU2 Project Website: https://su2code.github.io
+ *
+ * The SU2 Project is maintained by the SU2 Foundation
+ * (http://su2foundation.org)
+ *
+ * Copyright 2012-2020, 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 <http://www.gnu.org/licenses/>.
+ */
+
+#pragma once
+
+#include "CConductivityModel.hpp"
+
+/*!
+ * \class CConstantConductivity
+ * \brief Defines a constant thermal conductivity model.
+ * \author S.Vitale, M.Pini
+ */
+class CConstantConductivity final : public CConductivityModel {
+ public:
+  /*!
+   * \brief Constructor of the class.
+   */
+  CConstantConductivity(su2double kt_const) : kt_(kt_const) {}
+
+  /*!
+   * \brief return conductivity value.
+   */
+  su2double GetConductivity() const override { return kt_; }
+
+  /*!
+   * \brief return conductivity partial derivative value.
+   */
+  su2double Getdktdrho_T() const override { return dktdrho_t_; }
+
+  /*!
+   * \brief return conductivity partial derivative value.
+   */
+  su2double GetdktdT_rho() const override { return dktdt_rho_; }
+
+  /*!
+   * \brief Set thermal conductivity.
+   */
+  void SetConductivity(su2double t, su2double rho, su2double mu_lam, su2double mu_turb, su2double cp) override {}
+
+  /*!
+   * \brief Set thermal conductivity derivatives.
+   */
+  void SetDerConductivity(su2double t, su2double rho, su2double dmudrho_t, su2double dmudt_rho, su2double cp) override {
+  }
+
+ private:
+  su2double kt_{0.0};        /*!< \brief Thermal conductivity. */
+  su2double dktdrho_t_{0.0}; /*!< \brief DktDrho_T. */
+  su2double dktdt_rho_{0.0}; /*!< \brief DktDT_rho. */
+};