Files
CoolProp/dev/incompressible_liquids/CPIncomp/DataObjects.py
2019-01-12 20:48:56 -07:00

574 lines
23 KiB
Python

from __future__ import division, print_function
import numpy as np
import os, math
from .BaseObjects import IncompressibleData, IncompressibleFitter
from abc import ABCMeta
class SolutionData(object):
"""
A base class that defines all the variables needed
in order to make a proper fit. You can copy this code
put in your data and add some documentation for where the
information came from.
"""
ifrac_mass = "mass"
ifrac_mole = "mole"
ifrac_volume = "volume"
ifrac_undefined = "not defined"
ifrac_pure = "pure"
def __init__(self):
self.significantDigits = 7
self.name = None # Name of the current fluid
self.description = None # Description of the current fluid
self.reference = None # Reference data for the current fluid
self.Tmax = None # Maximum temperature in K
self.Tmin = None # Minimum temperature in K
self.xmax = 1.0 # Maximum concentration
self.xmin = 0.0 # Minimum concentration
self.xid = self.ifrac_undefined # Concentration is mole, mass or volume-based
self.TminPsat = None # Minimum saturation temperature in K
self.Tbase = None # Base value for temperature fits
self.xbase = None # Base value for concentration fits
self.temperature = IncompressibleData() # Temperature for data points in K
self.concentration = IncompressibleData() # Concentration data points in weight fraction
self.density = IncompressibleData() # Density in kg/m3
self.specific_heat = IncompressibleData() # Heat capacity in J/(kg.K)
self.viscosity = IncompressibleData() # Dynamic viscosity in Pa.s
self.conductivity = IncompressibleData() # Thermal conductivity in W/(m.K)
self.saturation_pressure = IncompressibleData() # Saturation pressure in Pa
self.T_freeze = IncompressibleData() # Freezing temperature in K
self.mass2input = IncompressibleData() # dd
self.volume2input = IncompressibleData() # dd
self.mole2input = IncompressibleData() # dd
# Some of the functions might need a guess array
#self.viscosity.type = self.viscosity.INCOMPRESSIBLE_EXPONENTIAL
#self.viscosity.coeffs = np.array([+7e+2, -6e+1, +1e+1])
#self.saturation_pressure.type = self.saturation_pressure.INCOMPRESSIBLE_EXPONENTIAL
#self.saturation_pressure.coeffs = np.array([-5e+3, +3e+1, -1e+1])
self.xref = 0.0
self.Tref = 0.0
self.pref = 0.0
self.href = 0.0
self.sref = 0.0
self.uref = 0.0
self.rhoref = 0.0
# def getDataObjects(self):
# objList = {}
# objList["temperature"] = self.temperature
# objList["concentration"] = self.concentration
# objList["density"] = self.density
# objList["specific_heat"] = self.specific_heat
# objList["viscosity"] = self.viscosity
# objList["conductivity"] = self.conductivity
# objList["saturation_pressure"] = self.saturation_pressure
# objList["T_freeze"] = self.T_freeze
# objList["volume2mass"] = self.volume2mass
# objList["mass2mole"] = self.mass2mole
# return objList
def roundSingle(self, x):
if x == 0.0: return 0.0
return round(x, self.significantDigits - int(math.floor(math.log10(abs(x)))) - 1)
def round(self, x):
r, c, res = IncompressibleFitter.shapeArray(x)
#digits = -1*np.floor(np.log10(res))+self.significantDigits-1
for i in range(r):
for j in range(c):
if np.isfinite(res[i, j]):
res[i, j] = self.roundSingle(res[i, j])
return res
# def getPolyObjects(self):
# objList = {}
# objList["density"] = self.density
# objList["specific heat"] = self.specific_heat
## objList["viscosity"] = self.viscosity
# objList["conductivity"] = self.conductivity
## objList["saturation_pressure"] = self.saturation_pressure
## objList["T_freeze"] = self.T_freeze
## objList["volume2mass"] = self.volume2mass
## objList["mass2mole"] = self.mass2mole
# return objList
#
# def getExpPolyObjects(self):
# objList = {}
# objList["viscosity"] = self.viscosity
# objList["saturation pressure"] = self.saturation_pressure
# return objList
def checkT(self, T, p, x):
if self.Tmin <= 0.: raise ValueError("Please specify the minimum temperature.")
if self.Tmax <= 0.: raise ValueError("Please specify the maximum temperature.");
if ((self.Tmin > T) or (T > self.Tmax)): raise ValueError("Your temperature {0} is not between {1} and {2}.".format(T, self.Tmin, self.Tmax))
TF = 0.0
if (self.T_freeze.type != IncompressibleData.INCOMPRESSIBLE_NOT_SET): TF = self.Tfreeze(T, p, x)
if (T < TF): raise ValueError("Your temperature {0} is below the freezing point of {1}.".format(T, TF))
else: return True
return False
def checkP(self, T, p, x):
ps = 0.0
if (self.saturation_pressure.type != IncompressibleData.INCOMPRESSIBLE_NOT_SET): ps = self.psat(T, p, x)
if (p < 0.0): raise ValueError("You cannot use negative pressures: {0} < {1}. ".format(p, 0.0))
if (p < ps): raise ValueError("Equations are valid for liquid phase only: {0} < {1}. ".format(p, ps))
else: return True
return False
def checkX(self, x):
if (self.xmin < 0.0 or self.xmin > 1.0): raise ValueError("Please specify the minimum concentration between 0 and 1.");
if (self.xmax < 0.0 or self.xmax > 1.0): raise ValueError("Please specify the maximum concentration between 0 and 1.");
if ((self.xmin > x) or (x > self.xmax)): raise ValueError("Your composition {0} is not between {1} and {2}.".format(x, self.xmin, self.xmax))
else: return True
return False
def checkTPX(self, T, p, x):
try:
return (self.checkT(T, p, x) and self.checkP(T, p, x) and self.checkX(x))
except ValueError as ve:
#print("Check failed: {0}".format(ve))
pass
return False
def rho(self, T, p=0.0, x=0.0, c=None):
if not self.checkTPX(T, p, x): return np.NAN
if c is None:
c = self.density.coeffs
if self.density.type == self.density.INCOMPRESSIBLE_POLYNOMIAL:
return np.polynomial.polynomial.polyval2d(T - self.Tbase, x - self.xbase, c)
else: raise ValueError("Unknown function.")
def c(self, T, p=0.0, x=0.0, c=None):
if not self.checkTPX(T, p, x): return np.NAN
if c is None:
c = self.specific_heat.coeffs
if self.specific_heat.type == self.specific_heat.INCOMPRESSIBLE_POLYNOMIAL:
return np.polynomial.polynomial.polyval2d(T - self.Tbase, x - self.xbase, c)
else: raise ValueError("Unknown function.")
def cp(self, T, p=0.0, x=0.0, c=None):
return self.c(T, p, x, c)
def cv(self, T, p=0.0, x=0.0, c=None):
return self.c(T, p, x, c)
def u(self, T, p=0.0, x=0.0, c=None):
if not self.checkTPX(T, p, x): return np.NAN
if c is None:
c = self.specific_heat.coeffs
if self.specific_heat.type == self.specific_heat.INCOMPRESSIBLE_POLYNOMIAL:
c_tmp = np.polynomial.polynomial.polyint(c)
return np.polynomial.polynomial.polyval2d(T - self.Tbase, x - self.xbase, c_tmp)
else: raise ValueError("Unknown function.")
# def u (T, p, x);
def h(self, T, p=0.0, x=0.0):
return self.h_u(T, p, x)
def visc(self, T, p=0.0, x=0.0, c=None):
if not self.checkTPX(T, p, x): return np.NAN
return self.viscosity.baseFunction(T, x, self.Tbase, self.xbase, c=c)
def cond(self, T, p=0.0, x=0.0, c=None):
if not self.checkTPX(T, p, x): return np.NAN
return self.conductivity.baseFunction(T, x, self.Tbase, self.xbase, c=c)
def psat(self, T, p=0.0, x=0.0, c=None):
if (T <= self.TminPsat): return 0.0
return self.saturation_pressure.baseFunction(T, x, self.Tbase, self.xbase, c=c)
def Tfreeze(self, T, p=0.0, x=0.0, c=None):
if c is None:
c = self.T_freeze.coeffs
if self.T_freeze.type == self.T_freeze.INCOMPRESSIBLE_POLYNOMIAL:
return np.polynomial.polynomial.polyval2d(p - 0.0, x - self.xbase, c)
elif self.T_freeze.type == self.T_freeze.INCOMPRESSIBLE_POLYOFFSET:
#if y!=0.0: raise ValueError("This is 1D only, use x not y.")
return self.T_freeze.basePolyOffset(c, x) # offset included in coeffs
elif self.T_freeze.type == self.T_freeze.INCOMPRESSIBLE_EXPONENTIAL:
#if y!=0.0: raise ValueError("This is 1D only, use x not y.")
return self.T_freeze.baseExponential(c, x)
elif self.T_freeze.type == IncompressibleData.INCOMPRESSIBLE_LOGEXPONENTIAL:
#if y!=0.0: raise ValueError("This is 1D only, use x not y.")
return self.T_freeze.baseLogexponential(c, x)
elif self.T_freeze.type == self.T_freeze.INCOMPRESSIBLE_EXPPOLYNOMIAL:
return np.exp(np.polynomial.polynomial.polyval2d(p - 0.0, x - self.xbase, c))
else:
raise ValueError("Unknown function: {0}.".format(self.T_freeze.type))
# def V2M (T, y);
# def M2M (T, x);
def h_u(self, T, p, x):
return self.u(T, p, x) + p / self.rho(T, p, x) - self.href
def u_h(self, T, p, x):
return self.h(T, p, x) - p / self.rho(T, p, x) + self.href
def set_reference_state(self, T0, p0, x0=0.0, h0=0.0, s0=0.0):
self.rhoref = self.rho(T0, p0, x0)
self.pref = p0
self.uref = h0 - p0 / self.rhoref
self.uref = self.u(T0, p0, x0)
self.href = h0
self.sref = s0
self.href = self.h(T0, p0, x0)
self.sref = self.s(T0, p0, x0)
class PureData(SolutionData):
"""
An extension of the solution data that makes it
easier to gather data for pure fluids.
"""
__metaclass__ = ABCMeta
def __init__(self):
SolutionData.__init__(self)
self.xbase = 0.0
self.xid = self.ifrac_pure
self.concentration.data = np.array([0]) # mass fraction
def reshapeData(self, dataArray, length):
"""
Makes any array 1-dimensional and implicitly
checks the length.
"""
if not dataArray is None:
return np.reshape(dataArray, (length, 1))
else:
return None
def reshapeAll(self):
len_T = len(self.temperature.data)
#len_x = len(self.concentration.data)
self.density.data = self.reshapeData(self.density.data, len_T)
self.specific_heat.data = self.reshapeData(self.specific_heat.data, len_T)
self.viscosity.data = self.reshapeData(self.viscosity.data, len_T)
self.conductivity.data = self.reshapeData(self.conductivity.data, len_T)
self.saturation_pressure.data = self.reshapeData(self.saturation_pressure.data, len_T)
#self.T_freeze.data = self.reshapeData(self.T_freeze.data, len_T)
#self.volume2mass.data = self.reshapeData(self.volume2mass.data, len_T)
#self.mass2mole.data = self.reshapeData(self.mass2mole.data, len_T)
class DigitalData(SolutionData):
"""
An extension of the solution data that makes it
easier to generate fitting data from fluids available
as Python packages.
"""
__metaclass__ = ABCMeta
def __init__(self):
SolutionData.__init__(self)
def getFile(self, data):
return os.path.join(os.path.dirname(__file__), 'data', self.name + "_" + data + ".txt")
def getFromFile(self, data):
fullPath = self.getFile(data)
_, _, res = IncompressibleFitter.shapeArray(np.loadtxt(fullPath))
return res
def writeToFile(self, data, array):
fullPath = self.getFile(data)
if not os.path.exists(os.path.dirname(fullPath)):
os.makedirs(os.path.dirname(fullPath))
stdFmt = "%1.{0}e".format(int(self.significantDigits - 1))
return np.savetxt(fullPath, array, fmt=stdFmt)
def getTrange(self):
if self.Tmin < self.Tmax:
return np.linspace(self.Tmin, self.Tmax, 20)
else:
return np.array([0.0])
def getxrange(self):
if self.xmin < self.xmax and self.xid != self.ifrac_pure:
return np.linspace(self.xmin, self.xmax, 20)
else:
return np.array([0.0])
def getArray(self, dataID=None, func=None, x_in=None, y_in=None, DEBUG=False):
""" Tries to read a data file, overwrites it if x or y do not match
:param dataID : ID to construct the path to the data file
:param func : Callable object that can take x_in and y_in
:param x_in : a 1D array in x direction or 2D with one column, most likely temperature
:param y_in : a 1D array in y direction or 2D with one row, most likely cocentration
:param DEBUG : a boolean that controls verbosity
:returns : Returns a tuple with three entries: x(1D),y(1D),data(2D)
"""
x = None
y = None
z = None
# First we try to read the file
if (not dataID is None and os.path.isfile(self.getFile(dataID))): # File found
fileArray = self.getFromFile(dataID)
x = np.copy(fileArray[1:, 0])
y = np.copy(fileArray[0, 1:])
z = np.copy(fileArray[1:, 1:])
else:
if DEBUG: print("No readable file found for {0}: {1}".format(dataID, self.getFile(dataID)))
updateFile = DEBUG
if not x_in is None: # Might need update
if not x is None: # Both given, check if different
mask = np.isfinite(x)
if IncompressibleFitter.allClose(x[mask], x_in[mask]):
if DEBUG: print("Both x-arrays are the same, no action required.")
updateFile = (updateFile or False) # Do not change a True value to False
else:
updateFile = True
if DEBUG: print("x-arrays do not match. {0} contains \n {1} \n and will be updated with \n {2}".format(self.getFile(dataID), x, x_in))
else: updateFile = True
elif x is None:
if DEBUG: print("Could not load x from file {0} and no x_in provided, aborting.".format(self.getFile(dataID)))
return None, None, None
else: updateFile = (updateFile or False) # Do not change a True value to False
if not y_in is None: # Might need update
if not y is None: # Both given, check if different
mask = np.isfinite(y)
if IncompressibleFitter.allClose(y[mask], y_in[mask]):
if DEBUG: print("Both y-arrays are the same, no action required.")
updateFile = (updateFile or False) # Do not change a True value to False
else:
updateFile = True
if DEBUG: print("y-arrays do not match. {0} contains \n {1} \n and will be updated with \n {2}".format(self.getFile(dataID), y, y_in))
else: updateFile = True
elif y is None:
if DEBUG: print("Could not load y from file {0} and no y_in provided, aborting.".format(self.getFile(dataID)))
return None, None, None
else: updateFile = (updateFile or False) # Do not change a True value to False
if DEBUG: print("Updating data file {0}".format(updateFile))
if not updateFile: return x, y, z # Done, data read from file
# Overwrite inputs
x = x_in
y = y_in
z = np.zeros((len(x) + 1, len(y) + 1))
r, c = z.shape
if func is None: raise ValueError("Need a function to update the data file.")
for i in range(r - 1):
for j in range(c - 1):
z[i + 1, j + 1] = func(x[i], y[j])
z[0, 0] = np.NaN
z[1:, 0] = x
z[0, 1:] = y
if not dataID is None:
self.writeToFile(dataID, z)
else:
if DEBUG: print("Not updating data file, dataID is missing.")
return x, y, z[1:, 1:]
class CoefficientData(SolutionData):
"""
A class to convert parameter arrays from different other sources
"""
__metaclass__ = ABCMeta
def __init__(self):
SolutionData.__init__(self)
self.reference = "Some other software"
def convertSecCoolArray(self, array):
if len(array) != 18:
raise ValueError("The length is not equal to 18!")
self.reference = "SecCool software"
array = np.array(array)
tmp = np.zeros((4, 6))
tmp[0][0] = array[0]
tmp[0][1] = array[1]
tmp[0][2] = array[2]
tmp[0][3] = array[3]
tmp[0][4] = array[4]
tmp[0][5] = array[5]
tmp[1][0] = array[6]
tmp[1][1] = array[7]
tmp[1][2] = array[8]
tmp[1][3] = array[9]
tmp[1][4] = array[10]
#tmp[1][5] = array[11]
tmp[2][0] = array[11]
tmp[2][1] = array[12]
tmp[2][2] = array[13]
tmp[2][3] = array[14]
#tmp[2][4] = array[4]
#tmp[2][5] = array[5]
tmp[3][0] = array[15]
tmp[3][1] = array[16]
tmp[3][2] = array[17]
#tmp[3][3] = array[3]
#tmp[3][4] = array[4]
#tmp[3][5] = array[5]
# Concentration is no longer handled in per cent!
for i in range(6):
tmp.T[i] *= 100.0**i
return tmp
def convertSecCoolTfreeze(self, array):
expo = -1.0
for i in range(len(array)):
array[i] = array[i] * np.power(100.0, expo + i)
array[1] = array[1] + 273.15
#self.T_freeze.type = self.T_freeze.INCOMPRESSIBLE_POLYOFFSET
return array
def convertMelinderArray(self, array):
"""The same function as the SecCool converter,
the original source code is slightly different though.
That is why the implementation is in a transposed form..."""
if len(array) != 18:
raise ValueError("The length is not equal to 18!")
#self.reference = "Melinder Book"
array = np.array(array)
tmp = np.zeros((6, 4))
tmp[0][0] = array[0]
tmp[0][1] = array[6]
tmp[0][2] = array[11]
tmp[0][3] = array[15]
tmp[1][0] = array[1]
tmp[1][1] = array[7]
tmp[1][2] = array[12]
tmp[1][3] = array[16]
tmp[2][0] = array[2]
tmp[2][1] = array[8]
tmp[2][2] = array[13]
tmp[2][3] = array[17]
tmp[3][0] = array[3]
tmp[3][1] = array[9]
tmp[3][2] = array[14]
tmp[4][0] = array[4]
tmp[4][1] = array[10]
tmp[5][0] = array[5]
# Concentration is no longer handled in per cent!
for i in range(6):
tmp[i] *= 100.0**i
return tmp.T
def convertMelinderMatrix(self, array):
"""Function to convert the full coefficient array
from the very first CoolProp implementation
based on the book by Melinder"""
if len(array) != 18:
raise ValueError("The length is not equal to 18!")
if len(array[0]) != 5:
raise ValueError("The length is not equal to 5!")
array = np.array(array)
tmp = np.zeros((18, 5))
for j in range(5):
tmp[0][j] = array[0][j]
tmp[1][j] = array[4][j]
tmp[2][j] = array[8][j]
tmp[3][j] = array[12][j]
tmp[4][j] = array[15][j]
tmp[5][j] = array[17][j]
tmp[6][j] = array[1][j]
tmp[7][j] = array[5][j]
tmp[8][j] = array[9][j]
tmp[9][j] = array[13][j]
tmp[10][j] = array[16][j]
tmp[11][j] = array[2][j]
tmp[12][j] = array[6][j]
tmp[13][j] = array[10][j]
tmp[14][j] = array[14][j]
tmp[15][j] = array[3][j]
tmp[16][j] = array[7][j]
tmp[17][j] = array[11][j]
return tmp
def setMelinderMatrix(self, matrix):
# matrix = np.array([
# [-26.29 , 958.1 ,3887 , 0.4175 , 1.153 ],
# [ -0.000002575 , -0.4151 , 7.201 , 0.0007271 , -0.03866 ],
# [ -0.000006732 , -0.002261 , -0.08979 , 0.0000002823 , 0.0002779 ],
# [ 0.000000163 , 0.0000002998 , -0.000439 , 0.000000009718 , -0.000001543 ],
# [ -1.187 , -1.391 , -18.5 , -0.004421 , 0.005448 ],
# [ -0.00001609 , -0.0151 , 0.2984 , -0.00002952 , 0.0001008 ],
# [ 0.000000342 , 0.0001113 , -0.001865 , 0.00000007336 , -0.000002809 ],
# [ 0.0000000005687, -0.0000003264 , -0.00001718 , 0.0000000004328 , 0.000000009811 ],
# [ -0.01218 , -0.01105 , -0.03769 , 0.00002044 , -0.0005552 ],
# [ 0.0000003865 , 0.0001828 , -0.01196 , 0.0000003413 , 0.000008384 ],
# [ 0.000000008768 , -0.000001641 , 0.00009801 , -0.000000003665 , -0.00000003997 ],
# [ -0.0000000002095, 0.0000000151 , 0.000000666 , -0.00000000002791 , -0.0000000003466 ],
# [ -0.00006823 , -0.0001208 , -0.003776 , 0.0000002943 , 0.000003038 ],
# [ 0.00000002137 , 0.000002992 , -0.00005611 , -0.0000000009646 , -0.00000007435 ],
# [ -0.0000000004271, 0.000000001455, -0.0000007811, 0.00000000003174 , 0.0000000007442 ],
# [ 0.0000001297 , 0.000004927 , -0.0001504 , -0.0000000008666 , 0.00000006669 ],
# [ -0.0000000005407, -0.0000001325 , 0.000007373 , -0.0000000000004573, -0.0000000009105 ],
# [ 0.00000002363 , -0.00000007727 , 0.000006433 , -0.0000000002033 , -0.0000000008472 ]
# ])
coeffs = self.convertMelinderMatrix(matrix).T
self.T_freeze.source = self.T_freeze.SOURCE_COEFFS
self.T_freeze.type = self.T_freeze.INCOMPRESSIBLE_POLYNOMIAL
self.T_freeze.coeffs = self.convertMelinderArray(coeffs[0])
self.T_freeze.coeffs[0, 0] += 273.15
self.T_freeze.coeffs = np.array([self.T_freeze.coeffs[0]])
# print(self.T_freeze.coeffs)
self.density.source = self.density.SOURCE_COEFFS
self.density.type = self.density.INCOMPRESSIBLE_POLYNOMIAL
self.density.coeffs = self.convertMelinderArray(coeffs[1])
self.specific_heat.source = self.specific_heat.SOURCE_COEFFS
self.specific_heat.type = self.specific_heat.INCOMPRESSIBLE_POLYNOMIAL
self.specific_heat.coeffs = self.convertMelinderArray(coeffs[2])
self.conductivity.source = self.conductivity.SOURCE_COEFFS
self.conductivity.type = self.conductivity.INCOMPRESSIBLE_POLYNOMIAL
self.conductivity.coeffs = self.convertMelinderArray(coeffs[3])
self.viscosity.source = self.viscosity.SOURCE_COEFFS
self.viscosity.type = self.viscosity.INCOMPRESSIBLE_EXPPOLYNOMIAL
self.viscosity.coeffs = self.convertMelinderArray(coeffs[4])
self.viscosity.coeffs[0, 0] -= math.log(1000) # Fixes the units mPa s -> Pa s