Ancillaries.h

Go to the documentation of this file.

#ifndef ANCILLARIES_H
#define ANCILLARIES_H
 
#include "CoolProp/Exceptions.h"
#include <utility>
#include <vector>
#include "CoolProp/detail/tools.h"  // for CoolPropDbl
#include "Eigen/Core"
#include "CoolProp/numerics/PolyMath.h"
 
namespace CoolProp {
 
class SurfaceTensionCorrelation
{
   public:
    std::vector<CoolPropDbl> a,  
      n,                         
      s;                         
    CoolPropDbl Tc;              
    std::size_t N;               
    std::string BibTeX;          
 
    SurfaceTensionCorrelation() : Tc(_HUGE), N(0) {}
 
    struct Values
    {
        std::vector<CoolPropDbl> a, n;
        CoolPropDbl Tc = _HUGE;
        std::string BibTeX;
    };
 
    explicit SurfaceTensionCorrelation(const Values& v)
      : a(v.a), n(v.n), s(v.n), Tc(v.Tc), N(v.n.size()), BibTeX(v.BibTeX) {}
 
    CoolPropDbl evaluate(CoolPropDbl T) {
        if (a.empty()) {
            throw NotImplementedError(format("surface tension curve not provided"));
        }
        if (T > Tc) {
            throw ValueError(format("Must be saturated state : T <= Tc"));
        }
        CoolPropDbl THETA = 1 - T / Tc;
        for (std::size_t i = 0; i < N; ++i) {
            s[i] = a[i] * pow(THETA, n[i]);
        }
        return std::accumulate(s.begin(), s.end(), 0.0);
    }
};
class SaturationAncillaryFunction
{
   public:
    enum ancillaryfunctiontypes
    {
        TYPE_NOT_SET = 0,
        TYPE_NOT_EXPONENTIAL,     
        TYPE_EXPONENTIAL,         
        TYPE_RATIONAL_POLYNOMIAL  
    };
 
    struct Values
    {
        ancillaryfunctiontypes type = TYPE_NOT_SET;
        Eigen::MatrixXd num_coeffs, den_coeffs;
        std::vector<double> n, t;
        CoolPropDbl max_abs_error = _HUGE;
        bool using_tau_r = false;
        CoolPropDbl reducing_value = _HUGE, T_r = _HUGE;
        CoolPropDbl Tmin = _HUGE, Tmax = _HUGE;
    };
 
   private:
    Eigen::MatrixXd num_coeffs,   
      den_coeffs;                 
    std::vector<double> n, t, s;  // For TYPE_NOT_EXPONENTIAL & TYPE_EXPONENTIAL
    union
    {
        CoolPropDbl max_abs_error;  
        struct
        {                                // For TYPE_NOT_EXPONENTIAL & TYPE_EXPONENTIAL
            bool using_tau_r;            
            CoolPropDbl reducing_value,  
              T_r;                       
            std::size_t N;               
        };
    };
    CoolPropDbl Tmax,  
      Tmin;            
    ancillaryfunctiontypes type;  
 
   public:
    SaturationAncillaryFunction()
      : type(TYPE_NOT_SET), Tmin(_HUGE), Tmax(_HUGE) {
 
        };
 
    explicit SaturationAncillaryFunction(const Values& v);
 
    bool enabled() {
        return type != TYPE_NOT_SET;
    }
 
    CoolPropDbl get_max_abs_error() {
        return max_abs_error;
    };
 
    double evaluate(double T);
 
    double invert(double value, double min_bound = -1, double max_bound = -1);
 
    double get_Tmin() {
        return Tmin;
    };
 
    double get_Tmax() {
        return Tmax;
    };
};
 
// ****************************************************************************
// ****************************************************************************
//                                 MELTING LINE
// ****************************************************************************
// ****************************************************************************
 
struct MeltingLinePiecewiseSimonSegment
{
    CoolPropDbl T_0, a, c, p_0, T_max, T_min, p_min, p_max;
};
struct MeltingLinePiecewiseSimonData
{
    std::vector<MeltingLinePiecewiseSimonSegment> parts;
};
 
class MeltingLinePiecewisePolynomialInTrSegment
{
   public:
    std::vector<CoolPropDbl> a, t;
    CoolPropDbl T_0, p_0, T_max, T_min, p_min, p_max;
    CoolPropDbl evaluate(CoolPropDbl T) {
        CoolPropDbl summer = 0;
        for (std::size_t i = 0; i < a.size(); ++i) {
            summer += a[i] * (pow(T / T_0, t[i]) - 1);
        }
        return p_0 * (1 + summer);
    }
};
struct MeltingLinePiecewisePolynomialInTrData
{
    std::vector<MeltingLinePiecewisePolynomialInTrSegment> parts;
};
 
class MeltingLinePiecewisePolynomialInThetaSegment
{
   public:
    std::vector<CoolPropDbl> a, t;
    CoolPropDbl T_0, p_0, T_max, T_min, p_min, p_max;
 
    CoolPropDbl evaluate(CoolPropDbl T) {
        CoolPropDbl summer = 0;
        for (std::size_t i = 0; i < a.size(); ++i) {
            summer += a[i] * pow(T / T_0 - 1, t[i]);
        }
        return p_0 * (1 + summer);
    }
};
struct MeltingLinePiecewisePolynomialInThetaData
{
    std::vector<MeltingLinePiecewisePolynomialInThetaSegment> parts;
};
 
class MeltingLineVariables
{
   public:
    enum MeltingLineVariablesEnum
    {
        MELTING_LINE_NOT_SET = 0,
        MELTING_LINE_SIMON_TYPE,                
        MELTING_LINE_POLYNOMIAL_IN_TR_TYPE,     
        MELTING_LINE_POLYNOMIAL_IN_THETA_TYPE,  
    };
    CoolPropDbl Tmin,  
      Tmax,            
      pmin,            
      pmax;            
 
    std::string BibTeX;                   
    CoolPropDbl T_m;                      
    MeltingLinePiecewiseSimonData simon;  
    MeltingLinePiecewisePolynomialInTrData
      polynomial_in_Tr;  
    MeltingLinePiecewisePolynomialInThetaData
      polynomial_in_Theta;  
    int type;
 
    MeltingLineVariables() : Tmin(_HUGE), Tmax(_HUGE), pmin(_HUGE), pmax(_HUGE), T_m(_HUGE), type(MELTING_LINE_NOT_SET) {};
 
    CoolPropDbl evaluate(int OF, int GIVEN, CoolPropDbl value);
 
    void set_limits();
 
    bool enabled() {
        return type != MELTING_LINE_NOT_SET;
    };
 
    [[nodiscard]] std::vector<std::pair<CoolPropDbl, CoolPropDbl>> get_parts_pranges() const;
};
 
} /* namespace CoolProp */
#endif