CoolProp fluid model

.github/workflows/regression.yml

@@ -26,7 +26,7 @@ jobs:
         config_set: [BaseMPI, ReverseMPI, ForwardMPI, BaseNoMPI, ReverseNoMPI, ForwardNoMPI, BaseOMP, ReverseOMP, ForwardOMP]
         include:
           - config_set: BaseMPI
-            flags: '-Denable-pywrapper=true -Denable-tests=true --warnlevel=3 --werror'
+            flags: '-Denable-pywrapper=true -Denable-coolprop=true -Denable-tests=true --warnlevel=2 --werror'
           - config_set: ReverseMPI
             flags: '-Denable-autodiff=true -Denable-normal=false -Denable-pywrapper=true -Denable-tests=true --warnlevel=3 --werror'
           - config_set: ForwardMPI
@@ -52,12 +52,12 @@ jobs:
           key: ${{ matrix.config_set }}-${{ github.sha }}
           restore-keys: ${{ matrix.config_set }}
       - name: Pre Cleanup 
-        uses: docker://ghcr.io/su2code/su2/build-su2:220614-1237
+        uses: docker://ghcr.io/su2code/su2/build-su2:221027-0442
         with:
           entrypoint: /bin/rm
           args: -rf install install_bin.tgz src ccache ${{ matrix.config_set }}
       - name: Build 
-        uses: docker://ghcr.io/su2code/su2/build-su2:220614-1237
+        uses: docker://ghcr.io/su2code/su2/build-su2:221027-0442
         with:
           args: -b ${{github.ref}} -f "${{matrix.flags}}"
       - name: Compress binaries
@@ -68,7 +68,7 @@ jobs:
           name: ${{ matrix.config_set }}
           path: install_bin.tgz
       - name: Post Cleanup
-        uses: docker://ghcr.io/su2code/su2/build-su2:220614-1237
+        uses: docker://ghcr.io/su2code/su2/build-su2:221027-0442
         with:
           entrypoint: /bin/rm
           args: -rf install install_bin.tgz src ccache ${{ matrix.config_set }}

.gitmodules

@@ -14,7 +14,10 @@
 	url = https://github.com/mesonbuild/meson
 [submodule "subprojects/Mutationpp"]
         path = subprojects/Mutationpp
-        url = https://github.com/mutationpp/Mutationpp.git	
+        url = https://github.com/mutationpp/Mutationpp.git
+[submodule "subprojects/CoolProp"]
+        path = subprojects/CoolProp
+        url = https://github.com/CoolProp/CoolProp.git
 [submodule "externals/opdi"]
 	path = externals/opdi
 	url = https://github.com/SciCompKL/OpDiLib

AUTHORS.md

@@ -109,6 +109,7 @@ Patrick Mischke
 Paul Urbanczyk
 Paul Zhang
 Pedro Gomes
+Peng Yan
 Pete Bachant
 RaulFeijo55
 Ruben Sanchez

Common/include/CConfig.hpp

@@ -80,6 +80,7 @@ class CConfig {
   su2double EA_ScaleFactor;       /*!< \brief Equivalent Area scaling factor */
   su2double AdjointLimit;         /*!< \brief Adjoint variable limit */
   string* ConvField;              /*!< \brief Field used for convergence check.*/
+  string FluidName;              /*!< \brief name of the applied fluid. */
 
   string* WndConvField;              /*!< \brief Function where to apply the windowed convergence criteria for the time average of the unsteady (single zone) flow problem. */
   unsigned short nConvField;         /*!< \brief Number of fields used to monitor convergence.*/
@@ -3752,6 +3753,11 @@ class CConfig {
    */
   unsigned short GetKind_FluidModel(void) const { return Kind_FluidModel; }
 
+  /*!
+   * \brief Returns the name of the fluid we are using in CoolProp.
+   */
+  string GetFluid_Name(void) const { return FluidName; }
+
   /*!
    * \brief Option to define the density model for incompressible flows.
    * \return Density model option

Common/include/option_structure.hpp

@@ -569,6 +569,7 @@ enum ENUM_FLUIDMODEL {
   MUTATIONPP = 7,         /*!< \brief Mutation++ gas model for nonequilibrium flow. */
   SU2_NONEQ = 8,          /*!< \brief User defined gas model for nonequilibrium flow. */
   FLUID_MIXTURE = 9,      /*!< \brief Species mixture model. */
+  COOLPROP = 10,          /*!< \brief Thermodynamics library. */
 };
 static const MapType<std::string, ENUM_FLUIDMODEL> FluidModel_Map = {
   MakePair("STANDARD_AIR", STANDARD_AIR)
@@ -581,6 +582,7 @@ static const MapType<std::string, ENUM_FLUIDMODEL> FluidModel_Map = {
   MakePair("MUTATIONPP", MUTATIONPP)
   MakePair("SU2_NONEQ", SU2_NONEQ)
   MakePair("FLUID_MIXTURE", FLUID_MIXTURE)
+  MakePair("COOLPROP", COOLPROP)
 };
 
 /*!

Common/src/CConfig.cpp

@@ -1150,6 +1150,8 @@ void CConfig::SetConfig_Options() {
   /*!\par CONFIG_CATEGORY: FluidModel \ingroup Config*/
   /*!\brief FLUID_MODEL \n DESCRIPTION: Fluid model \n OPTIONS: See \link FluidModel_Map \endlink \n DEFAULT: STANDARD_AIR \ingroup Config*/
   addEnumOption("FLUID_MODEL", Kind_FluidModel, FluidModel_Map, STANDARD_AIR);
+  /*!\brief FLUID_NAME \n DESCRIPTION: Fluid name \n OPTIONS: see coolprop homepage \n DEFAULT: nitrogen \ingroup Config*/
+  addStringOption("FLUID_NAME", FluidName, string("nitrogen"));
 
 
   /*!\par CONFIG_CATEGORY: Freestream Conditions \ingroup Config*/
@@ -3383,6 +3385,11 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i
   }
 #endif
 
+  /*--- Check if CoolProp is used with non-dimensionalization. ---*/
+  if (Kind_FluidModel == COOLPROP && Ref_NonDim != DIMENSIONAL) {
+    SU2_MPI::Error("CoolProp can not be used with non-dimensionalization.", CURRENT_FUNCTION);
+  }
+
   /*--- STL_BINARY output not implemented yet, but already a value in option_structure.hpp---*/
   for (unsigned short iVolumeFile = 0; iVolumeFile < nVolumeOutputFiles; iVolumeFile++) {
     if (VolumeOutputFiles[iVolumeFile] == OUTPUT_TYPE::STL_BINARY){

SU2_CFD/include/fluid/CCoolProp.hpp

@@ -0,0 +1,148 @@
+/*!
+ * \file CCoolProp.hpp
+ * \brief Defines the state-of-the-art fluid model from CoolProp library.
+ * \author P. Yan, G. Gori, A. Guardone
+ * \version 7.4.0 "Blackbird"
+ *
+ * SU2 Project Website: https://su2code.github.io
+ *
+ * The SU2 Project is maintained by the SU2 Foundation
+ * (http://su2foundation.org)
+ *
+ * Copyright 2012-2022, 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 "CFluidModel.hpp"
+#if defined(HAVE_COOLPROP) && !defined(CODI_FORWARD_TYPE) && !defined(CODI_REVERSE_TYPE)
+#define USE_COOLPROP
+namespace CoolProp {
+  class AbstractState;
+}
+#endif
+#include <memory>
+
+
+/*!
+ * \class CCoolProp
+ * \brief Child class for defining fluid model from CoolProp library.
+ * \author: P.Yan
+ */
+class CCoolProp final : public CFluidModel {
+ private:
+  su2double Gamma{1.4};           /*!< \brief Ratio of Specific Heats. */
+  su2double Gas_Constant{297};    /*!< \brief specific Gas Constant. */
+  su2double Pressure_Critical{0.0};   /*!< \brief critical pressure */
+  su2double Temperature_Critical{0.0};    /*!< \brief critical temperature */
+  su2double acentric_factor{0.0};      /*!< \brief acentric factor */
+#ifdef USE_COOLPROP
+  std::unique_ptr<CoolProp::AbstractState> fluid_entity;   /*!< \brief fluid entity */
+#endif
+
+ public:
+  /*!
+   * \brief Constructor of the class.
+   */
+  CCoolProp(string fluidname);
+
+#ifdef USE_COOLPROP
+  /*!
+   * \brief Destructor of the class.
+   * \note Needs to be defined in the .cpp to allow using only a forward declaration of CoolProp::AbstractState.
+   */
+  ~CCoolProp();
+
+  /*!
+   * \brief Set the Dimensionless State using Density and Internal Energy
+   * \param[in] rho - first thermodynamic variable.
+   * \param[in] e - second thermodynamic variable.
+   */
+  void SetTDState_rhoe(su2double rho, su2double e) override;
+
+  /*!
+   * \brief Set the Dimensionless State using Pressure  and Temperature
+   * \param[in] P - first thermodynamic variable.
+   * \param[in] T - second thermodynamic variable.
+   */
+  void SetTDState_PT(su2double P, su2double T) override;
+
+  /*!
+   * \brief Set the Dimensionless State using Pressure and Density
+   * \param[in] P - first thermodynamic variable.
+   * \param[in] rho - second thermodynamic variable.
+   */
+  void SetTDState_Prho(su2double P, su2double rho) override;
+
+  /*!
+   * \brief Set the Dimensionless Internal Energy using Pressure and Density
+   * \param[in] P - first thermodynamic variable.
+   * \param[in] rho - second thermodynamic variable.
+   */
+  void SetEnergy_Prho(su2double P, su2double rho) override;
+
+  /*!
+   * \brief Set the Dimensionless State using Enthalpy and Entropy
+   * \param[in] th1 - first thermodynamic variable (h).
+   * \param[in] th2 - second thermodynamic variable (s).
+   */
+  void SetTDState_hs(su2double h, su2double s) override;
+
+  /*!
+   * \brief Set the Dimensionless State using Density and Temperature
+   * \param[in] th1 - first thermodynamic variable (rho).
+   * \param[in] th2 - second thermodynamic variable (T).
+   */
+  void SetTDState_rhoT(su2double rho, su2double T) override;
+
+  /*!
+   * \brief Set the Dimensionless State using Pressure and Entropy
+   * \param[in] th1 - first thermodynamic variable (P).
+   * \param[in] th2 - second thermodynamic variable (s).
+   */
+  void SetTDState_Ps(su2double P, su2double s) override;
+
+  /*!
+   * \brief compute some derivatives of enthalpy and entropy needed for subsonic inflow BC
+   * \param[in] th1 - first thermodynamic variable (P).
+   * \param[in] th2 - second thermodynamic variable (rho).
+   */
+  void ComputeDerivativeNRBC_Prho(su2double P, su2double rho) override;
+#endif
+
+  /*!
+   * \brief Get the value of the critical pressure.
+   * \return Critical pressure.
+   */
+  su2double GetPressure_Critical(void) const { return Pressure_Critical; }
+
+  /*!
+   * \brief Get the value of the critical temperature.
+   * \return Critical temperature.
+   */
+  su2double GetTemperature_Critical(void) const { return Temperature_Critical; }
+
+  /*!
+   * \brief Get the value of specific gas constant.
+   * \return Value of the constant: Gamma
+   */
+  su2double GetGas_Constant(void) const { return Gas_Constant; }
+
+  /*!
+   * \brief Get the value of specific gas constant.
+   * \return Value of the constant: Gamma
+   */
+  su2double GetGamma(void) const { return Gamma; }
+};

SU2_CFD/obj/Makefile.am

@@ -54,6 +54,7 @@ libSU2Core_sources = ../src/definition_structure.cpp \
   ../src/fluid/CFluidScalar.cpp \
   ../src/fluid/CPengRobinson.cpp \
   ../src/fluid/CVanDerWaalsGas.cpp \
+  ../src/fluid/CCoolProp.cpp \
   ../src/fluid/CNEMOGas.cpp \
   ../src/fluid/CSU2TCLib.cpp \
   ../src/fluid/CMutationTCLib.cpp \
@@ -246,6 +247,11 @@ su2_cfd_cxx_flags += @MUTATIONPP_CXX@
 su2_cfd_ldadd += @MUTATIONPP_LD@
 # endif
 
+
+su2_cfd_cxx_flags += @CoolProp_CXX@
+su2_cfd_ldadd += @CoolProp_LD@
+
+
 if BUILD_NORMAL
 libSU2Core_a_SOURCES = $(libSU2Core_sources)
 libSU2Core_a_CXXFLAGS = $(libSU2Core_cxx_flags)

SU2_CFD/src/fluid/CCoolProp.cpp

@@ -0,0 +1,119 @@
+/*!
+ * \file CCoolProp.cpp
+ * \brief Source of the fluid model from CoolProp.
+ * \author P. Yan, G. Gori, A. Guardone
+ * \version 7.4.0 "Blackbird"
+ *
+ * SU2 Project Website: https://su2code.github.io
+ *
+ * The SU2 Project is maintained by the SU2 Foundation
+ * (http://su2foundation.org)
+ *
+ * Copyright 2012-2022, 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/>.
+ */
+
+#include "../../include/fluid/CCoolProp.hpp"
+
+#ifdef USE_COOLPROP
+#include "CoolProp.h"
+#include "AbstractState.h"
+
+CCoolProp::CCoolProp(string fluidname) : CFluidModel() {
+  fluid_entity = std::unique_ptr<CoolProp::AbstractState>(CoolProp::AbstractState::factory("HEOS",fluidname));
+  Gas_Constant = fluid_entity->gas_constant()/fluid_entity->molar_mass();
+  Pressure_Critical = fluid_entity->p_critical();
+  Temperature_Critical = fluid_entity->T_critical();
+  acentric_factor = fluid_entity->acentric_factor();
+}
+
+CCoolProp::~CCoolProp() {}
+
+void CCoolProp::SetTDState_rhoe(su2double rho, su2double e) {
+  //cout<<"p "<<Pressure<<"Pc "<<Pressure_Critical<<"T "<<Temperature<<"Tc"<<Temperature_Critical<<endl;
+  Density = rho;
+  StaticEnergy = e;
+  fluid_entity->update(CoolProp::DmassUmass_INPUTS, Density, StaticEnergy);
+  Cp = fluid_entity->cpmass();
+  Cv = fluid_entity->cvmass();
+  Gamma = Cp / Cv;
+  Pressure = fluid_entity->p();
+  Temperature = fluid_entity->T();
+  Entropy = fluid_entity->smass();
+  dPdrho_e = fluid_entity->first_partial_deriv(CoolProp::iP, CoolProp::iDmass, CoolProp::iUmass);
+  dPde_rho = fluid_entity->first_partial_deriv(CoolProp::iP, CoolProp::iUmass, CoolProp::iDmass);
+  dTdrho_e = fluid_entity->first_partial_deriv(CoolProp::iT, CoolProp::iDmass, CoolProp::iUmass);
+  dTde_rho = fluid_entity->first_partial_deriv(CoolProp::iT, CoolProp::iUmass, CoolProp::iDmass);
+  if (fluid_entity->phase() == 6) {
+      fluid_entity->specify_phase(CoolProp::iphase_gas);
+      SetTDState_PT(Pressure,Temperature);
+  }
+  else{
+      SoundSpeed2 = pow(fluid_entity->speed_sound(), 2);
+  }
+}
+
+void CCoolProp::SetTDState_PT(su2double P, su2double T) {
+  fluid_entity->update(CoolProp::PT_INPUTS, P, T);
+  su2double rho = fluid_entity->rhomass();
+  su2double e = fluid_entity->umass();
+  SetTDState_rhoe(rho, e);
+}
+
+void CCoolProp::SetTDState_Prho(su2double P, su2double rho) {
+  fluid_entity->update(CoolProp::DmassP_INPUTS, rho, P);
+  su2double e = fluid_entity->umass();
+  SetTDState_rhoe(rho, e);
+}
+
+void CCoolProp::SetEnergy_Prho(su2double P, su2double rho) {
+  fluid_entity->update(CoolProp::DmassP_INPUTS, rho, P);
+  StaticEnergy = fluid_entity->umass();
+}
+
+void CCoolProp::SetTDState_hs(su2double h, su2double s) {
+  fluid_entity->update(CoolProp::HmassSmass_INPUTS, h, s);
+  su2double rho = fluid_entity->rhomass();
+  su2double e = fluid_entity->umass();
+  SetTDState_rhoe(rho, e);
+}
+
+void  CCoolProp::SetTDState_Ps(su2double P, su2double s) {
+  fluid_entity->update(CoolProp::PSmass_INPUTS, P, s);
+  su2double Rho = fluid_entity->rhomass();
+  su2double e = fluid_entity->umass();
+  SetTDState_rhoe(Rho, e);
+}
+
+void CCoolProp::SetTDState_rhoT(su2double rho, su2double T) {
+  fluid_entity->update(CoolProp::DmassT_INPUTS, rho, T);
+  su2double Rho = fluid_entity->rhomass();
+  su2double e = fluid_entity->umass();
+  SetTDState_rhoe(Rho, e);
+}
+
+void  CCoolProp::ComputeDerivativeNRBC_Prho(su2double P, su2double rho) {
+  SetTDState_Prho(P, rho);
+  dhdrho_P = fluid_entity->first_partial_deriv(CoolProp::iHmass,CoolProp::iDmass,CoolProp::iP);
+  dhdP_rho = fluid_entity->first_partial_deriv(CoolProp::iHmass,CoolProp::iP,CoolProp::iDmass);
+  dsdP_rho = fluid_entity->first_partial_deriv(CoolProp::iSmass,CoolProp::iP,CoolProp::iDmass);
+  dsdrho_P = fluid_entity->first_partial_deriv(CoolProp::iSmass,CoolProp::iDmass,CoolProp::iP);
+}
+
+#else
+CCoolProp::CCoolProp(string fluidname) {
+  SU2_MPI::Error("SU2 was not compiled with CoolProp (-Denable-coolprop=true). Note that CoolProp cannot be used with directdiff or autodiff", CURRENT_FUNCTION);
+}
+#endif