From 00e9ddf29165a6d5488ae06ca7db390c13b47511 Mon Sep 17 00:00:00 2001 From: jowr Date: Sun, 29 Jun 2014 18:03:37 +0200 Subject: [PATCH] Added 2D function to JSON readers, fixed 2D polynomial fitting --- dev/IncompressibleLiquids/LinearAlgebra.ipynb | 306 +++++++++++++ dev/IncompressibleLiquids/PureExample.json | 57 +++ dev/IncompressibleLiquids/SecCoolExample.json | 182 ++++++++ .../SolutionExample.json | 63 +++ .../data_incompressible.py | 418 +++++++++++++++++- .../fit_incompressible.py | 122 +++++ .../fit_incompressible_simple.py | 333 -------------- .../json_incompressible.py | 104 +++++ externals/eigen | 2 +- include/IncompressibleFluid.h | 4 +- include/rapidjson/rapidjson_include.h | 68 +++ .../Incompressible/IncompressibleLibrary.cpp | 39 +- 12 files changed, 1324 insertions(+), 374 deletions(-) create mode 100644 dev/IncompressibleLiquids/LinearAlgebra.ipynb create mode 100644 dev/IncompressibleLiquids/PureExample.json create mode 100644 dev/IncompressibleLiquids/SecCoolExample.json create mode 100644 dev/IncompressibleLiquids/SolutionExample.json delete mode 100644 dev/IncompressibleLiquids/fit_incompressible_simple.py create mode 100644 dev/IncompressibleLiquids/json_incompressible.py diff --git a/dev/IncompressibleLiquids/LinearAlgebra.ipynb b/dev/IncompressibleLiquids/LinearAlgebra.ipynb new file mode 100644 index 00000000..6c511307 --- /dev/null +++ b/dev/IncompressibleLiquids/LinearAlgebra.ipynb @@ -0,0 +1,306 @@ +{ + "metadata": { + "name": "", + "signature": "sha256:6de5b361bf724471b38e3a32d973641aeac571541b15819486c1f0449e310b07" + }, + "nbformat": 3, + "nbformat_minor": 0, + "worksheets": [ + { + "cells": [ + { + "cell_type": "heading", + "level": 1, + "metadata": {}, + "source": [ + "Linear algebra for multidimensional polynomial fitting" + ] + }, + { + "cell_type": "code", + "collapsed": false, + "input": [ + "import numpy as np\n", + "import itertools\n", + "# A small data set\n", + "#T_in = np.array([20,24,35,40])+273.15\n", + "#x_in = np.array([47,55,70,78,82])/100.0\n", + "#rho_in = np.array([[1047,1033,1020,1000,990],[997,983,970,950,940],[947,933,920,900,895],[987,883,870,850,845]])\n", + "#\n", + "# large data set\n", + "T_in = np.array([ -45 , -40 , -35 , -30 , -25 , -20 , -15 , -10])+273.15 # Kelvin\n", + "x_in = np.array([ 5 , 10 , 15 , 20 , 25 , 30 , 35 ])/100.0 # mass fraction\n", + " \n", + "rho_in = np.array([\n", + " [1064.0, 1054.6, 1045.3, 1036.3, 1027.4, 1018.6, 1010.0],\n", + " [1061.3, 1052.1, 1043.1, 1034.3, 1025.6, 1017.0, 1008.6],\n", + " [1057.6, 1048.8, 1040.1, 1031.5, 1023.1, 1014.8, 1006.7],\n", + " [1053.1, 1044.6, 1036.2, 1028.0, 1019.9, 1012.0, 1004.1],\n", + " [1047.5, 1039.4, 1031.5, 1023.7, 1016.0, 1008.4, 1000.9],\n", + " [1040.7, 1033.2, 1025.7, 1018.4, 1011.2, 1004.0, 997.0],\n", + " [1032.3, 1025.3, 1018.5, 1011.7, 1005.1, 998.5, 992.0],\n", + " [1021.5, 1015.3, 1009.2, 1003.1, 997.1, 991.2, 985.4]]) # kg/m3" + ], + "language": "python", + "metadata": {}, + "outputs": [], + "prompt_number": 13 + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Minimizing the squared error $\\epsilon(\\mathbf{c}) = \\sqrt{\\sum (\\mathbf{z} - \\mathbf{A} \\cdot \\mathbf{c})^2 }$ can be achieved by solving the system of orthogonal equations given by $\\mathbf{A}^\\text{T}\\mathbf{A} \\cdot \\mathbf{c} =\\mathbf{A}^\\text{T}\\mathbf{z}$. Using Python tools, we can leave this to the software and we only have to construct the Van der Monde matrix $\\mathbf{A}$ of the independent variable and equate it with the result vector of the dependent variable." + ] + }, + { + "cell_type": "code", + "collapsed": false, + "input": [ + "def getCoeffs1d(x,z,order):\n", + " if (len(x)" + ] + } + ], + "prompt_number": 47 + }, + { + "cell_type": "code", + "collapsed": false, + "input": [], + "language": "python", + "metadata": {}, + "outputs": [] + } + ], + "metadata": {} + } + ] +} \ No newline at end of file diff --git a/dev/IncompressibleLiquids/PureExample.json b/dev/IncompressibleLiquids/PureExample.json new file mode 100644 index 00000000..72ff7e62 --- /dev/null +++ b/dev/IncompressibleLiquids/PureExample.json @@ -0,0 +1,57 @@ +{ + "T_freeze": { + "coeffs": "null", + "type": 0 + }, + "Tbase": 0.0, + "Tmax": 533.15, + "Tmin": 188.14999999999998, + "TminPsat": 313.15, + "conductivity": { + "coeffs": "null", + "type": 0 + }, + "density": { + "coeffs": [ + [ + -5.439005439004754e-06 + ], + [ + 0.005552564102563395 + ], + [ + -2.6422007556329925 + ], + [ + 1197.8367716325445 + ] + ], + "type": "polynomial" + }, + "description": "Heat transfer fluid TherminolD12 by Solutia", + "mass2mole": { + "coeffs": "null", + "type": 0 + }, + "name": "PureExample", + "reference": "Solutia data sheet", + "saturation_pressure": { + "coeffs": "null", + "type": 0 + }, + "specific_heat": { + "coeffs": "null", + "type": 0 + }, + "viscosity": { + "coeffs": "null", + "type": 0 + }, + "volume2mass": { + "coeffs": "null", + "type": 0 + }, + "xbase": 0.0, + "xmax": null, + "xmin": null +} \ No newline at end of file diff --git a/dev/IncompressibleLiquids/SecCoolExample.json b/dev/IncompressibleLiquids/SecCoolExample.json new file mode 100644 index 00000000..4bec2ca3 --- /dev/null +++ b/dev/IncompressibleLiquids/SecCoolExample.json @@ -0,0 +1,182 @@ +{ + "T_freeze": { + "coeffs": [ + 27.7555556, + -22.9732217, + -1.10405072, + -0.0120762281, + -9.343458e-05 + ], + "type": "polynomial" + }, + "Tbase": 268.66999999999996, + "Tmax": 293.15, + "Tmin": 223.14999999999998, + "TminPsat": 293.15, + "conductivity": { + "coeffs": [ + [ + 0.40820667, + -0.003981687, + 1.583368e-05, + -3.552049e-07, + -9.884176e-10, + 4.46e-10 + ], + [ + 0.0006629321, + -2.686475e-05, + 9.03915e-07, + -2.128257e-08, + -5.562e-10, + 0.0 + ], + [ + 3.685975e-07, + 7.188416e-08, + -1.041773e-08, + 2.278001e-10, + 0.0, + 0.0 + ], + [ + 4.703395e-08, + 7.612361e-11, + -2.734e-10, + 0.0, + 0.0, + 0.0 + ] + ], + "type": "polynomial" + }, + "density": { + "coeffs": [ + [ + 960.246658, + 1.29038391, + -0.016104252, + -0.0001969888, + 1.131559e-05, + 9.181999e-08 + ], + [ + -0.402034827, + -0.0162463989, + 0.0001623301, + 4.367343e-06, + 1.199e-08, + 0.0 + ], + [ + -0.0025204776, + 0.0001101514, + -2.320217e-07, + 7.794999e-08, + 0.0, + 0.0 + ], + [ + 9.937483e-06, + -1.346886e-06, + 4.141999e-08, + 0.0, + 0.0, + 0.0 + ] + ], + "type": "polynomial" + }, + "description": "Methanol solution", + "mass2mole": { + "coeffs": "null", + "type": 0 + }, + "name": "SecCoolExample", + "reference": "SecCool software", + "saturation_pressure": { + "coeffs": "null", + "type": 0 + }, + "specific_heat": { + "coeffs": [ + [ + 3822.97123, + -23.1224095, + 0.0678775826, + 0.0022413893, + -0.0003045332, + -4.758e-06 + ], + [ + 2.35014495, + 0.178883941, + 0.0006828, + 0.0002101166, + -9.812e-06, + 0.0 + ], + [ + -0.0004724176, + -0.0003317949, + 0.0001002032, + -5.306e-06, + 0.0, + 0.0 + ], + [ + 4.242194e-05, + 2.34719e-05, + -1.894e-06, + 0.0, + 0.0, + 0.0 + ] + ], + "type": "polynomial" + }, + "viscosity": { + "coeffs": [ + [ + 1.47255255, + 0.0022218998, + -0.0004406139, + 6.047984e-06, + -1.95473e-07, + -2.372e-09 + ], + [ + -0.0411841566, + 0.0001784479, + -3.564413e-06, + 4.064671e-08, + 1.915e-08, + 0.0 + ], + [ + 0.0002572862, + -9.226343e-07, + -2.178577e-08, + -9.529999e-10, + 0.0, + 0.0 + ], + [ + -1.699844e-06, + -1.023552e-07, + 4.482e-09, + 0.0, + 0.0, + 0.0 + ] + ], + "type": "polynomial" + }, + "volume2mass": { + "coeffs": "null", + "type": 0 + }, + "xbase": 0.3157, + "xmax": 0.5, + "xmin": 0.0 +} \ No newline at end of file diff --git a/dev/IncompressibleLiquids/SolutionExample.json b/dev/IncompressibleLiquids/SolutionExample.json new file mode 100644 index 00000000..e77242ac --- /dev/null +++ b/dev/IncompressibleLiquids/SolutionExample.json @@ -0,0 +1,63 @@ +{ + "T_freeze": { + "coeffs": "null", + "type": 0 + }, + "Tbase": 0.0, + "Tmax": null, + "Tmin": null, + "TminPsat": null, + "conductivity": { + "coeffs": "null", + "type": 0 + }, + "density": { + "coeffs": [ + [ + -1.1876190476190525, + -1.0978571428571455, + -1.008571428571432, + -0.926190476190484, + -0.8433333333333355, + -0.7626190476190485, + -0.6845238095238149 + ], + [ + 1338.9886190476204, + 1308.851107142858, + 1278.9555714285727, + 1250.893690476193, + 1222.8398333333341, + 1195.3998690476192, + 1168.7407738095253 + ] + ], + "type": "polynomial" + }, + "description": "Ethanol ice slurry", + "mass2mole": { + "coeffs": "null", + "type": 0 + }, + "name": "SolutionExample", + "reference": "SecCool software", + "saturation_pressure": { + "coeffs": "null", + "type": 0 + }, + "specific_heat": { + "coeffs": "null", + "type": 0 + }, + "viscosity": { + "coeffs": "null", + "type": 0 + }, + "volume2mass": { + "coeffs": "null", + "type": 0 + }, + "xbase": 0.0, + "xmax": null, + "xmin": null +} \ No newline at end of file diff --git a/dev/IncompressibleLiquids/data_incompressible.py b/dev/IncompressibleLiquids/data_incompressible.py index 18128de1..61eb4414 100644 --- a/dev/IncompressibleLiquids/data_incompressible.py +++ b/dev/IncompressibleLiquids/data_incompressible.py @@ -1,4 +1,4 @@ -import numpy +import numpy,itertools,scipy.interpolate class LiquidData(object): @@ -532,29 +532,405 @@ class ZS55(LiquidData): #dataObj = Therminol72() #for i in range(len(dataObj.T)): # print str(dataObj.T[i])+", "+str(numpy.log(dataObj.mu_dyn[i])) + +class IncompressibleData(object): + def __init__(self): + self.INCOMPRESSIBLE_NOT_SET = 0 + self.INCOMPRESSIBLE_POLYNOMIAL = 'polynomial' + self.INCOMPRESSIBLE_EXPONENTIAL = 'exponential' + self.INCOMPRESSIBLE_EXPPOLYNOMIAL = 'exppolynomial' + self.type = self.INCOMPRESSIBLE_NOT_SET + self.coeffs = None #numpy.zeros((4,4)) + self.data = None #numpy.zeros((10,10)) + + def fit(self, T, x=0.0, Tbase=0.0, xbase=0.0): + (cr,cc) = self.coeffs.shape + (dr,dc) = self.data.shape + (Tr,Tc) = (len(T),1) #T.shape #(len(T),1) + if (Tc!=1): raise ValueError("Temperature has to be a 2D array with one column.") + if (Tr!=dr): raise ValueError("Temperature and fitting data have to have the same number of rows.") + + if self.type==self.INCOMPRESSIBLE_POLYNOMIAL: + if (numpy.max(x)==0.0): # x not given + x = numpy.array([0.0]) + if (cc==1): # requires 1D coefficients + if (dc==1): + #z = numpy.polyfit(T, self.data[:,0], cr-1) + #self.coeffs = numpy.zeros((cr,1)) + #self.coeffs[:,0] = z[::-1] + #print self.coeffs + self.coeffs = self.getCoeffs2d(T-Tbase, x-xbase, self.data, cr-1, cc-1) + #print self.coeffs + else: + raise ValueError("Cannot use 2D data with 1D references") + else: + raise ValueError("Cannot use 2D coefficients without concentration input") + else: # Assume 2D input + (xr,xc) = (1,len(x))#x.shape#(1,len(x)) + if (xr!=1): raise ValueError("Concentration has to be a 2D array with one row.") + if (xc!=dc): raise ValueError("Concentration and fitting data have to have the same number of columns.") + #if (cr!=cc): raise ValueError("For now, you have to have the same number of exponents for x and y.") + self.coeffs = self.getCoeffs2d(T-Tbase, x-xbase, self.data, cr-1, cc-1) + #raise ValueError("Cannot use 2D coefficients yet") + print self.coeffs + else: + raise ValueError("Cannot fit that function.") + -class SolutionData(LiquidData): + def getCoeffs1d(self, x,z,order): + if (len(x) > coeffs; diff --git a/include/rapidjson/rapidjson_include.h b/include/rapidjson/rapidjson_include.h index 00bb6b81..a4f0434f 100644 --- a/include/rapidjson/rapidjson_include.h +++ b/include/rapidjson/rapidjson_include.h @@ -100,6 +100,46 @@ namespace cpjson return out; }; + /// A convenience function to get a 2D double array compactly + inline std::vector< std::vector > get_double_array2D(rapidjson::Value &v) + { + std::vector< std::vector > out; + std::vector tmp; + if (!v.IsArray()) { throw CoolProp::ValueError("input is not an array"); } + for (rapidjson::Value::ValueIterator itr = v.Begin(); itr != v.End(); ++itr) + { + if (!itr->IsArray()) { throw CoolProp::ValueError("input is not a 2D array"); } + tmp.clear(); + for (rapidjson::Value::ValueIterator i = itr->Begin(); i != itr->End(); ++i) + { + if (!i->IsNumber()){throw CoolProp::ValueError("input is not a number");} + tmp.push_back(i->GetDouble()); + } + out.push_back(tmp); + } + return out; + }; + + /// A convenience function to get a 2D long double array compactly + inline std::vector< std::vector > get_long_double_array2D(rapidjson::Value &v) + { + std::vector< std::vector > out; + std::vector tmp; + if (!v.IsArray()) { throw CoolProp::ValueError("input is not an array"); } + for (rapidjson::Value::ValueIterator itr = v.Begin(); itr != v.End(); ++itr) + { + if (!itr->IsArray()) { throw CoolProp::ValueError("input is not a 2D array"); } + tmp.clear(); + for (rapidjson::Value::ValueIterator i = itr->Begin(); i != itr->End(); ++i) + { + if (!i->IsNumber()){throw CoolProp::ValueError("input is not a number");} + tmp.push_back(i->GetDouble()); + } + out.push_back(tmp); + } + return out; + }; + /// A convenience function to get a long double array compactly inline std::vector get_long_double_array(rapidjson::Value &v, std::string name) { @@ -135,6 +175,34 @@ namespace cpjson return buffer.GetString(); }; +// /// A convenience function to set an array compactly +// template +// inline void set_array(const char *key, const std::vector &vec, rapidjson::Value &value, rapidjson::Document &doc) +// { +// rapidjson::Value _v(rapidjson::kArrayType); +// for (unsigned int i = 0; i < vec.size(); ++i) +// { +// _v.PushBack(vec[i],doc.GetAllocator()); +// } +// value.AddMember(key, _v, doc.GetAllocator()); +// }; +// +// /// A convenience function to set an array compactly +// template +// inline void set_array2D(const char *key, const std::vector< std::vector > &vec, rapidjson::Value &value, rapidjson::Document &doc) +// { +// rapidjson::Value _i(rapidjson::kArrayType); +// rapidjson::Value _j; +// for (unsigned int i = 0; i < vec.size(); ++i) { +// _j = rapidjson::Value(rapidjson::kArrayType); +// for (unsigned int j = 0; j < vec[i].size(); ++j) { +// _j.PushBack(vec[j],doc.GetAllocator()); +// } +// _i.PushBack(_j,doc.GetAllocator()); +// } +// value.AddMember(key, _i, doc.GetAllocator()); +// }; + /// A convenience function to set a double array compactly inline void set_double_array(const char *key, const std::vector &vec, rapidjson::Value &value, rapidjson::Document &doc) { diff --git a/src/Backends/Incompressible/IncompressibleLibrary.cpp b/src/Backends/Incompressible/IncompressibleLibrary.cpp index 0bc9ebcf..528e9f79 100644 --- a/src/Backends/Incompressible/IncompressibleLibrary.cpp +++ b/src/Backends/Incompressible/IncompressibleLibrary.cpp @@ -11,28 +11,33 @@ IncompressibleData JSONIncompressibleLibrary::parse_coefficients(rapidjson::Valu if (obj.HasMember(id.c_str())) { //rapidjson::Value value = obj[id.c_str()]; if (obj[id.c_str()].HasMember("type")){ - std::string type = cpjson::get_string(obj[id.c_str()], "type"); - if (!type.compare("polynomial")){ - fluidData.type = CoolProp::IncompressibleData::INCOMPRESSIBLE_POLYNOMIAL; - fluidData.coeffs = vec_to_eigen(cpjson::get_double_array(obj[id.c_str()]["coeffs"])); - return fluidData; - } - else if (!type.compare("exponential")){ - fluidData.type = CoolProp::IncompressibleData::INCOMPRESSIBLE_EXPONENTIAL; - fluidData.coeffs = vec_to_eigen(cpjson::get_double_array(obj[id.c_str()]["coeffs"])); - return fluidData; - } - else if (!type.compare("exppolynomial")){ - fluidData.type = CoolProp::IncompressibleData::INCOMPRESSIBLE_EXPPOLYNOMIAL; - fluidData.coeffs = vec_to_eigen(cpjson::get_double_array(obj[id.c_str()]["coeffs"])); - return fluidData; + if (obj[id.c_str()].HasMember("coeffs")){ + std::string type = cpjson::get_string(obj[id.c_str()], "type"); + if (!type.compare("polynomial")){ + fluidData.type = CoolProp::IncompressibleData::INCOMPRESSIBLE_POLYNOMIAL; + fluidData.coeffs = vec_to_eigen(cpjson::get_double_array2D(obj[id.c_str()]["coeffs"])); + return fluidData; + } + else if (!type.compare("exponential")){ + fluidData.type = CoolProp::IncompressibleData::INCOMPRESSIBLE_EXPONENTIAL; + fluidData.coeffs = vec_to_eigen(cpjson::get_double_array(obj[id.c_str()]["coeffs"])); + return fluidData; + } + else if (!type.compare("exppolynomial")){ + fluidData.type = CoolProp::IncompressibleData::INCOMPRESSIBLE_EXPPOLYNOMIAL; + fluidData.coeffs = vec_to_eigen(cpjson::get_double_array2D(obj[id.c_str()]["coeffs"])); + return fluidData; + } + else{ + throw ValueError(format("The type [%s] is not understood for [%s] of incompressible fluids. Please check your JSON file.", type.c_str(), id.c_str())); + } } else{ - throw ValueError(format("The type [%s] is not understood for [%s] of incompressible fluids. Please check your JSON file.", type.c_str(), id.c_str())); + throw ValueError(format("Your file does not have an entry for \"coeffs\" in [%s], which is vital for this function.", id.c_str())); } } else{ - throw ValueError(format("Your file does not have an entry for \"type\" of [%s], which is vital for this function.", id.c_str())); + throw ValueError(format("Your file does not have an entry for \"type\" in [%s], which is vital for this function.", id.c_str())); } } else{