Merge branch 'master' into draw_process

This commit is contained in:
Jorrit Wronski
2015-10-08 13:33:25 +02:00
20 changed files with 832 additions and 97 deletions

View File

@@ -1270,27 +1270,32 @@ if (COOLPROP_SNIPPETS)
foreach (snippet ${snippets})
get_filename_component(snippet_name ${snippet} NAME)
get_filename_component(snippet_exe ${snippet} NAME_WE)
message(STATUS "snippet_name = ${snippet_name}")
add_executable (${snippet_name} ${snippet})
add_dependencies (${snippet_name} CoolProp)
target_link_libraries (${snippet_name} CoolProp)
add_executable (${snippet_exe} ${snippet})
add_dependencies (${snippet_exe} CoolProp)
target_link_libraries (${snippet_exe} CoolProp)
if(UNIX)
target_link_libraries (${snippet_name} ${CMAKE_DL_LIBS})
target_link_libraries (${snippet_exe} ${CMAKE_DL_LIBS})
endif()
if ( MSVC )
set_target_properties( ${snippet_name} PROPERTIES RUNTIME_OUTPUT_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/bin )
set_target_properties( ${snippet_name} PROPERTIES RUNTIME_OUTPUT_DIRECTORY_DEBUG ${CMAKE_CURRENT_BINARY_DIR}/bin )
set_target_properties( ${snippet_name} PROPERTIES RUNTIME_OUTPUT_DIRECTORY_RELEASE ${CMAKE_CURRENT_BINARY_DIR}/bin )
set_target_properties( ${snippet_exe} PROPERTIES RUNTIME_OUTPUT_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/bin )
set_target_properties( ${snippet_exe} PROPERTIES RUNTIME_OUTPUT_DIRECTORY_DEBUG ${CMAKE_CURRENT_BINARY_DIR}/bin )
set_target_properties( ${snippet_exe} PROPERTIES RUNTIME_OUTPUT_DIRECTORY_RELEASE ${CMAKE_CURRENT_BINARY_DIR}/bin )
# etc for the other available configuration types (MinSizeRel, RelWithDebInfo)
set(BIN_PATH "${CMAKE_CURRENT_BINARY_DIR}/bin")
else()
set(BIN_PATH "${CMAKE_CURRENT_BINARY_DIR}")
endif ()
# Run it and save the output to a file with .output appended
message(STATUS "${CMAKE_CURRENT_BINARY_DIR}/bin/${snippet_name} > ${CMAKE_CURRENT_BINARY_DIR}/bin/${snippet_name}.output")
SET_PROPERTY(TARGET ${snippet_exe} APPEND_STRING PROPERTY COMPILE_FLAGS " -DEXTERNC")
add_custom_command(TARGET ${snippet_name}
# Run it and save the output to a file with .output appended
add_custom_command(TARGET ${snippet_exe}
POST_BUILD
COMMAND ${CMAKE_CURRENT_BINARY_DIR}/bin/${snippet_name} > ${CMAKE_CURRENT_SOURCE_DIR}/Web/coolprop/snippets/${snippet_name}.output)
COMMAND ${BIN_PATH}/${snippet_exe} > ${CMAKE_CURRENT_SOURCE_DIR}/Web/coolprop/snippets/${snippet_name}.output)
endforeach()
@@ -1302,7 +1307,7 @@ if (COOLPROP_CLANG_ADDRESS_SANITIZER)
list(APPEND APP_SOURCES "${CMAKE_CURRENT_SOURCE_DIR}/src/Tests/test_main.cxx")
# CATCH TEST, compile everything with catch and set test entry point
add_executable (CatchTestRunner ${APP_SOURCES})
add_dependencies (CatchTestRuanner generate_headers)
add_dependencies (CatchTestRunner generate_headers)
set_target_properties (CatchTestRunner PROPERTIES COMPILE_FLAGS "${COMPILE_FLAGS} -DENABLE_CATCH")
set(CMAKE_EXE_LINKER_FLAGS "-fsanitize=address -lstdc++")
if(UNIX)

View File

@@ -216,7 +216,7 @@ More Information
The tables are stored in a zipped format using the msgpack package and miniz. If you want to see what data is serialized in the tabular data, you can unzip and unpack into python (or other high-level languages) using something roughly like::
import msgpack, zlib, StringIO
import msgpack, zlib, StringIO, numpy as np
with open(r'/path/to/home/.CoolProp/Tables/HelmholtzEOSBackend(R245fa)/single_phase_logph.bin.z','rb') as fp:
ph = zlib.decompress(fp.read())

View File

@@ -1,7 +1 @@
Loading table: C:\Users\Belli/.CoolProp/Tables/HelmholtzEOSBackend(Water)/single_phase_logph.bin.z
Loaded table: C:\Users\Belli/.CoolProp/Tables/HelmholtzEOSBackend(Water)/single_phase_logph.bin.z in 0.139 sec.
Loading table: C:\Users\Belli/.CoolProp/Tables/HelmholtzEOSBackend(Water)/single_phase_logpT.bin.z
Loaded table: C:\Users\Belli/.CoolProp/Tables/HelmholtzEOSBackend(Water)/single_phase_logpT.bin.z in 0.133 sec.
Loading table: C:\Users\Belli/.CoolProp/Tables/HelmholtzEOSBackend(Water)/pure_saturation.bin.z
Loaded table: C:\Users\Belli/.CoolProp/Tables/HelmholtzEOSBackend(Water)/pure_saturation.bin.z in 0.003 sec.
value(all): 8339004.420912, 0.45 us/call
value(all): 8339004.432517, 0.59425 us/call

View File

@@ -1,5 +1,5 @@
#include "CoolProp.h"
#include "MixtureDerivatives.h"
#include "Backends/Helmholtz/MixtureDerivatives.h"
#include <iostream>
using namespace CoolProp;
int main()
@@ -11,7 +11,7 @@ int main()
shared_ptr<HelmholtzEOSMixtureBackend> HEOS(new HelmholtzEOSMixtureBackend(components));
HelmholtzEOSMixtureBackend &rHEOS = *(HEOS.get());
HEOS->set_mole_fractions(z);
HEOS->specify_phase(iphase_gas); // So that we don't do a
HEOS->specify_phase(iphase_gas); // So that we don't do a phase check
HEOS->update(DmolarT_INPUTS, 300, 300);
std::vector<std::string> terms;
@@ -43,7 +43,7 @@ int main()
for (std::vector<std::string>::iterator it = terms.begin(); it != terms.end(); ++it)
{
if (!it->compare("p")){
printf("p: %0.16Lg\n", HEOS->p());
printf("p: %0.16g\n", HEOS->p());
}
else if (!it->compare("p2(deriv)")){
printf("p calculated by rho*R*T*(1+delta*deltadar_dDelta): %0.16Lg\n", HEOS->rhomolar()*HEOS->gas_constant()*HEOS->T()*(1+HEOS->delta()*HEOS->dalphar_dDelta()));
@@ -52,16 +52,16 @@ int main()
printf("dalphar_dDelta: %0.16Lg\n", HEOS->dalphar_dDelta());
}
else if (!it->compare("rhor")){
printf("rhor: %0.16Lg\n", HEOS->get_reducing_state().rhomolar);
printf("rhor: %0.16g\n", HEOS->get_reducing_state().rhomolar);
}
else if (!it->compare("Tr")){
printf("Tr: %0.16Lg\n", HEOS->get_reducing_state().T);
printf("Tr: %0.16g\n", HEOS->get_reducing_state().T);
}
else if (!it->compare("dTr_dxi")){
printf("dTr_dxi: %0.16Lg\n", HEOS->Reducing.p->dTrdxi__constxj(rHEOS.get_mole_fractions(), 0, XN_DEPENDENT));
printf("dTr_dxi: %0.16Lg\n", HEOS->Reducing->dTrdxi__constxj(rHEOS.get_mole_fractions(), 0, XN_DEPENDENT));
}
else if (!it->compare("drhor_dxi")){
printf("drhor_dxi: %0.16Lg\n", HEOS->Reducing.p->drhormolardxi__constxj(rHEOS.get_mole_fractions(), 0, XN_DEPENDENT));
printf("drhor_dxi: %0.16Lg\n", HEOS->Reducing->drhormolardxi__constxj(rHEOS.get_mole_fractions(), 0, XN_DEPENDENT));
}
else if(!it->compare("ndpdV__constT_n")){
printf("ndpdV__constT_n: %0.16Lg\n", MixtureDerivatives::ndpdV__constT_n(rHEOS));

View File

@@ -1,25 +1,25 @@
p: 675616.5260238834
p calculated by rho*R*T*(1+delta*deltadar_dDelta): 675616.5260238834
p: 675615.7215706832
p calculated by rho*R*T*(1+delta*deltadar_dDelta): 675615.7215706832
rhor: 5375.073740431392
Tr: 354.7802371525603
dalphar_dDelta: -1.740349537597104
dTr_dxi: -60.84527088778825
drhor_dxi: 1587.916430678275
ndpdV__constT_n: -181257322.2040499
dpdxj__constT_V_xi(0): 48181.54565089135
dalphar_dxi|T,V,xk(0): 0.06467412878843236
dalphar_dxi|tau_delta_xk(0): 4.187330988393492e-005
ln_fugacity_coefficient(0): -0.04441530290159221
dTr_dxi: -60.84527088778821
drhor_dxi: 1587.916430678276
ndpdV__constT_n: -181257106.3818679
dpdxj__constT_V_xi(0): 48181.48828136717
dalphar_dxi|T,V,xk(0): 0.06467412878843232
dalphar_dxi|tau_delta_xk(0): 4.187330988391106e-05
ln_fugacity_coefficient(0): -0.04441530290159228
ln_fugacity_coefficient(1): -0.1090894316900246
ndpdni__constT_V_nj(0): 640327.2332516684
tau*d_ndalphardni_dTau(0): -0.09709340835471267
tau*d_ndalphardni_dTau(1): -0.2484618994956092
delta*d_ndalphardni_dDelta(0): -0.04715912844791718
ndpdni__constT_V_nj(0): 640326.4708172516
tau*d_ndalphardni_dTau(0): -0.09709340835471265
tau*d_ndalphardni_dTau(1): -0.2484618994956091
delta*d_ndalphardni_dDelta(0): -0.04715912844791722
delta*d_ndalphardni_dDelta(1): -0.1115469220723501
d_ndalphardni_dxj__constdelta_tau_xi(0, 0): -0.04493881035511202
d_ndalphardni_dxj__constT_V_xi(0, 0): -0.01435530858047472
d_ndalphardni_dxj__constdelta_tau_xi(0, 0): -0.04493881035511207
d_ndalphardni_dxj__constT_V_xi(0, 0): -0.01435530858047477
dln_fugacity_coefficient_dxj__constT_p_xi(0,0): -0.01791995312054628
dln_fugacity_coefficient_dxj__constT_p_xi(1,0): 0.005973317706848738
d2nalphar_dxj_dni__constT_V(0,0): 0.05031882020795764
dln_fugacity_coefficient_dxj__constT_p_xi(1,0): 0.005973317706848762
d2nalphar_dxj_dni__constT_V(0,0): 0.05031882020795755
d2nalphar_dxj_dni__constT_V(1,0): 0
delta*d2alphar_dxi_dDelta(0): 0.001221607211340736
delta*d2alphar_dxi_dDelta(0): 0.00122160721133991

View File

@@ -12,9 +12,10 @@ int main()
// Second type (C++ only, a bit faster, allows for vector inputs and outputs)
std::vector<std::string> fluids; fluids.push_back("Propane"); fluids.push_back("Ethane");
std::vector<std::string> outputs; outputs.push_back("Dmolar");
std::vector<double> T(1,298), p(1e5);
std::cout << PropsSImulti(outputs,"T", T, "P", p, "", fluids, z)[0][0] << std::endl;
std::vector<double> T(1,298), p(1,1e5);
std::cout << PropsSImulti(outputs,"T", T, "P", p, "", fluids, z)[0][0] << std::endl; // Default backend is HEOS
std::cout << PropsSImulti(outputs,"T", T, "P", p, "HEOS", fluids, z)[0][0] << std::endl;
// Comment me out if REFPROP is not installed
std::cout << PropsSImulti(outputs,"T", T, "P", p, "REFPROP", fluids, z)[0][0] << std::endl;
return EXIT_SUCCESS;

View File

@@ -63,6 +63,13 @@ Octave Requirements
* SWIG
* C#
OSX
---
For OSX, to install the necessary tools using homebrew, you can do::
homebrew install mono
Linux
-----

View File

@@ -33,4 +33,23 @@ You can build the static library using::
# Build the makefile using CMake
cmake .. -DCOOLPROP_STATIC_LIBRARY=ON
# Make the static library
make VERBOSE=1
cmake --build .
Usage
-----
On linux and OSX, you can use the compiled static library in your application by doing something like this (starting in the directory ``build`` relative to the root of the source as in the above compilation step)::
g++ -lCoolProp -ldl -L. -I../../include main.cpp
where main.cpp could have the contents for instance of::
#include "CoolProp.h"
#include <iostream>
int main()
{
std::cout << CoolProp::PropsSI("T","P",101325,"Q",0,"Water") << std::endl;
return 1;
}

View File

@@ -0,0 +1,491 @@
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Helmholtz energy conversion of cubics"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"As in Michelsen's book, a generalized volume-translated cubic EOS can be given the common form:\n",
"$$ p = \\frac{RT}{v-b}-\\frac{a(T)}{(v+\\Delta_1b)(v+\\Delta_2b)} $$\n",
"\n",
"$$ p = \\frac{RT}{v+c-b}-\\frac{a(T)}{(v+c+\\Delta_1b)(v+c+\\Delta_2b)} $$\n"
]
},
{
"cell_type": "code",
"execution_count": 78,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"IPython console for SymPy 0.7.6 (Python 2.7.10-32-bit) (ground types: python)\n"
]
}
],
"source": [
"import sys\n",
"from __future__ import division\n",
"from sympy import *\n",
"from IPython.display import display, Math, Latex\n",
"from IPython.core.display import display_html\n",
"init_session(quiet=True)\n",
"init_printing()"
]
},
{
"cell_type": "code",
"execution_count": 79,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"b,rho_r,rho,tau,delta,T_r,R,Delta_1,Delta_2,T,v,c = symbols('b,rho_r,rho,tau,delta,T_r,R,Delta_1,Delta_2,T,v,c', positive = True)\n",
"a = symbols('a', cls=Function, positive = True)(tau)"
]
},
{
"cell_type": "code",
"execution_count": 95,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"normal_cubic = True\n",
"cubic = 'VdW'"
]
},
{
"cell_type": "code",
"execution_count": 96,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/latex": [
"$$p = \\frac{R T}{- b + c + v} - \\frac{a{\\left (\\tau \\right )}}{\\left(c + v\\right)^{2}}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"execution_count": 96,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"p = R*T/(v+c-b)-a/((v+c+Delta_1*b)*(v+c+Delta_2*b))\n",
"if cubic == 'SRK':\n",
" p = p.subs(Delta_1,1).subs(Delta_2,0)\n",
"elif cubic == 'VdW':\n",
" p = p.subs(Delta_1,0).subs(Delta_2,0)\n",
"elif cubic == 'PR':\n",
" p = p.subs(Delta_1,1+sqrt(2)).subs(Delta_2,1-sqrt(2))\n",
"Math('p = '+latex(p))"
]
},
{
"cell_type": "code",
"execution_count": 97,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/latex": [
"$$Z = \\frac{v}{R T} \\left(\\frac{R T}{- b + c + v} - \\frac{a{\\left (\\tau \\right )}}{\\left(c + v\\right)^{2}}\\right)$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"execution_count": 97,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Z = p*v/(R*T)\n",
"Math('Z = '+latex(Z))"
]
},
{
"cell_type": "code",
"execution_count": 98,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/latex": [
"$$Z = \\frac{1}{R T \\rho} \\left(\\frac{R T}{- b + c + \\frac{1}{\\rho}} - \\frac{a{\\left (\\tau \\right )}}{\\left(c + \\frac{1}{\\rho}\\right)^{2}}\\right)$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$Z = \\frac{1}{R T \\delta \\rho_{r}} \\left(\\frac{R T}{- b + c + \\frac{1}{\\delta \\rho_{r}}} - \\frac{a{\\left (\\tau \\right )}}{\\left(c + \\frac{1}{\\delta \\rho_{r}}\\right)^{2}}\\right)$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{d\\alpha^r}{d\\delta}|\\tau = - \\frac{a{\\left (\\tau \\right )}}{R T c^{2} \\delta^{2} \\rho_{r} + 2 R T c \\delta + \\frac{R T}{\\rho_{r}}} + \\frac{1}{- b \\delta^{2} \\rho_{r} + c \\delta^{2} \\rho_{r} + \\delta} - \\frac{1}{\\delta}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"execution_count": 98,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Z = Z.replace(v, 1/rho)\n",
"display(Math('Z = '+latex(Z)))\n",
"#Z = Z.replace(T, T_r/tau)\n",
"Z = Z.replace(rho, delta*rho_r)\n",
"display(Math('Z = '+latex(Z)))\n",
"dalphar_dDelta = (Z-1)/delta\n",
"dalphar_dDelta = expand(simplify(dalphar_dDelta))\n",
"Math('\\\\frac{d\\\\alpha^r}{d\\delta}|\\\\tau = ' + latex(dalphar_dDelta))"
]
},
{
"cell_type": "code",
"execution_count": 99,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/latex": [
"$$\\alpha^r = \\frac{1}{R T c \\left(c \\delta \\rho_{r} + 1\\right)} \\left(- R T c \\left(c \\delta \\rho_{r} + 1\\right) \\log{\\left (\\frac{\\delta \\rho_{r} \\left(b - c\\right) - 1}{\\rho_{r} \\left(b - c\\right)} \\right )} + \\left(c \\delta \\rho_{r} + 1\\right) \\left(R T c \\log{\\left (- \\frac{1}{\\rho_{r} \\left(b - c\\right)} \\right )} - a{\\left (\\tau \\right )}\\right) + a{\\left (\\tau \\right )}\\right)$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\lim_{c \\to 0}\\alpha^r = - \\log{\\left (b \\delta \\rho_{r} - 1 \\right )} + i \\pi - \\frac{\\delta \\rho_{r}}{R T} a{\\left (\\tau \\right )}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"antideriv_pieces = 0\n",
"for arg in expand(simplify(dalphar_dDelta))._args:\n",
" antideriv_pieces += integrate(arg,delta)\n",
"\n",
"antideriv_pieces = simplify(antideriv_pieces)\n",
"alphar = antideriv_pieces - antideriv_pieces.subs(delta,0)\n",
"\n",
"display(Math('\\\\alpha^r = ' + latex(simplify(alphar))))\n",
"display(Math('\\\\lim_{c \\\\to 0}\\\\alpha^r = ' + latex(simplify(limit(alphar,c,0)))))"
]
},
{
"cell_type": "code",
"execution_count": 100,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# Maybe turn off volume translation\n",
"if normal_cubic:\n",
" alphar = limit(alphar,c,0)"
]
},
{
"cell_type": "code",
"execution_count": 101,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"def format_deriv(arg, itau, idel, RHS):\n",
" \"\"\" \n",
" A function for giving a nice latex representation of \n",
" the partial derivative in question \n",
" \"\"\"\n",
" if itau+idel == 1:\n",
" numexp = ''\n",
" else:\n",
" numexp = '^{{{s:d}}}'.format(s=itau+idel)\n",
" \n",
" if itau == 0:\n",
" tau = ''\n",
" elif itau == 1:\n",
" tau = '\\\\partial \\\\tau'\n",
" else:\n",
" tau = '\\\\partial \\\\tau^{{{s:d}}}'.format(s=itau)\n",
" \n",
" if idel == 0:\n",
" delta = ''\n",
" elif idel == 1:\n",
" delta = '\\\\partial \\\\delta'\n",
" else:\n",
" delta = '\\\\partial \\\\delta^{{{s:d}}}'.format(s=idel)\n",
" \n",
" temp = '\\\\frac{{\\\\partial{{{numexp:s}}} {arg:s}}}{{{{{tau:s}}}{{{delta:s}}}}} = '\n",
" return Math(temp.format(**locals()) + latex(RHS))"
]
},
{
"cell_type": "code",
"execution_count": 102,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/latex": [
"$$\\frac{\\partial{} \\alpha^r}{{}{\\partial \\delta}} = \\frac{1}{R T} \\left(- \\frac{R T}{\\delta - \\frac{1}{b \\rho_{r}}} - \\rho_{r} a{\\left (\\tau \\right )}\\right)$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{} \\alpha^r}{{\\partial \\tau}{}} = - \\frac{\\delta \\rho_{r}}{R T} \\frac{d}{d \\tau} a{\\left (\\tau \\right )}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{2}} \\alpha^r}{{}{\\partial \\delta^{2}}} = \\frac{1}{\\left(\\delta - \\frac{1}{b \\rho_{r}}\\right)^{2}}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{2}} \\alpha^r}{{\\partial \\tau}{\\partial \\delta}} = - \\frac{\\rho_{r}}{R T} \\frac{d}{d \\tau} a{\\left (\\tau \\right )}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{2}} \\alpha^r}{{\\partial \\tau^{2}}{}} = - \\frac{\\delta \\rho_{r}}{R T} \\frac{d^{2}}{d \\tau^{2}} a{\\left (\\tau \\right )}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{3}} \\alpha^r}{{}{\\partial \\delta^{3}}} = - \\frac{2}{\\left(\\delta - \\frac{1}{b \\rho_{r}}\\right)^{3}}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{3}} \\alpha^r}{{\\partial \\tau}{\\partial \\delta^{2}}} = 0$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{3}} \\alpha^r}{{\\partial \\tau^{2}}{\\partial \\delta}} = - \\frac{\\rho_{r}}{R T} \\frac{d^{2}}{d \\tau^{2}} a{\\left (\\tau \\right )}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{3}} \\alpha^r}{{\\partial \\tau^{3}}{}} = - \\frac{\\delta \\rho_{r}}{R T} \\frac{d^{3}}{d \\tau^{3}} a{\\left (\\tau \\right )}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{4}} \\alpha^r}{{}{\\partial \\delta^{4}}} = \\frac{6}{\\left(\\delta - \\frac{1}{b \\rho_{r}}\\right)^{4}}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{4}} \\alpha^r}{{\\partial \\tau}{\\partial \\delta^{3}}} = 0$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{4}} \\alpha^r}{{\\partial \\tau^{2}}{\\partial \\delta^{2}}} = 0$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{4}} \\alpha^r}{{\\partial \\tau^{3}}{\\partial \\delta}} = - \\frac{\\rho_{r}}{R T} \\frac{d^{3}}{d \\tau^{3}} a{\\left (\\tau \\right )}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"$$\\frac{\\partial{^{4}} \\alpha^r}{{\\partial \\tau^{4}}{}} = - \\frac{\\delta \\rho_{r}}{R T} \\frac{d^{4}}{d \\tau^{4}} a{\\left (\\tau \\right )}$$"
],
"text/plain": [
"<IPython.core.display.Math object>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"for deriv_count in range(1,5):\n",
" for dtau in range(deriv_count+1):\n",
" ddelta = deriv_count-dtau\n",
" #print dtau, ddelta\n",
" display(format_deriv('\\\\alpha^r', dtau, ddelta, diff(diff(alphar,tau,dtau),delta,ddelta)))"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 2",
"language": "python",
"name": "python2"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.10"
}
},
"nbformat": 4,
"nbformat_minor": 0
}

View File

@@ -414,9 +414,23 @@ public:
/// Get the mole fractions of the equilibrium liquid phase
std::vector<CoolPropDbl> mole_fractions_liquid(void){ return calc_mole_fractions_liquid(); };
/// Get the mole fractions of the equilibrium liquid phase (but as a double for use in SWIG wrapper)
std::vector<double> mole_fractions_liquid_double(void){
std::vector<CoolPropDbl> x = calc_mole_fractions_liquid();
return std::vector<double>(x.begin(), x.end());
};
/// Get the mole fractions of the equilibrium vapor phase
std::vector<CoolPropDbl> mole_fractions_vapor(void){ return calc_mole_fractions_vapor(); };
/// Get the mole fractions of the equilibrium vapor phase (but as a double for use in SWIG wrapper)
std::vector<double> mole_fractions_vapor_double(void){
std::vector<CoolPropDbl> y = calc_mole_fractions_vapor();
return std::vector<double>(y.begin(), y.end());
};
/// Get the mole fractions of the components
virtual const std::vector<CoolPropDbl> & get_mole_fractions(void) = 0;
/// Update the state using two state variables
virtual void update(CoolProp::input_pairs input_pair, double Value1, double Value2) = 0;
@@ -424,9 +438,6 @@ public:
/// Some or all of the guesses will be used - this is backend dependent
virtual void update_with_guesses(CoolProp::input_pairs input_pair, double Value1, double Value2, const GuessesStructure &guesses){ throw NotImplementedError("update_with_guesses is not implemented for this backend"); };
/// Get the mole fractions of the components
virtual const std::vector<CoolPropDbl> & get_mole_fractions(void) = 0;
/// A function that says whether the backend instance can be instantiated in the high-level interface
/// In general this should be true, except for some other backends (especially the tabular backends)
/// To disable use in high-level interface, implement this function and return false

View File

@@ -74,7 +74,7 @@ AbstractState * AbstractState::factory(const std::string &backend, const std::ve
{
throw ValueError("TREND backend not yet implemented");
}
else if (!backend.compare("?"))
else if (!backend.compare("?") || backend.empty())
{
std::size_t idel = fluid_names[0].find("::");
// Backend has not been specified, and we have to figure out what the backend is by parsing the string

View File

@@ -28,7 +28,7 @@ IncompressibleBackend::IncompressibleBackend(IncompressibleFluid* fluid) {
this->fluid = fluid;
this->set_reference_state();
if (this->fluid->is_pure()){
this->set_fractions(std::vector<CoolPropDbl>(1,0.0));
this->set_fractions(std::vector<CoolPropDbl>(1,1.0));
}
}
@@ -36,7 +36,7 @@ IncompressibleBackend::IncompressibleBackend(const std::string &fluid_name) {
this->fluid = &get_incompressible_fluid(fluid_name);
this->set_reference_state();
if (this->fluid->is_pure()){
this->set_fractions(std::vector<CoolPropDbl>(1,0.0));
this->set_fractions(std::vector<CoolPropDbl>(1,1.0));
}
}
@@ -135,6 +135,14 @@ bool IncompressibleBackend::clear() {
this->_rhomass.clear();
this->_smass.clear();
this->_umass.clear();
this->_drhodTatPx.clear();
this->_dsdTatPx.clear();
this->_dhdTatPx.clear();
this->_dsdTatPxdT.clear();
this->_dhdTatPxdT.clear();
this->_dsdpatTx.clear();
this->_dhdpatTx.clear();
// Done
return true;
}
@@ -267,6 +275,35 @@ double IncompressibleBackend::cmass(void){
return _cmass;
}
double IncompressibleBackend::drhodTatPx(void){
if (!_drhodTatPx) _drhodTatPx = calc_drhodTatPx(_T,_p,_fractions[0]);
return _drhodTatPx;
}
double IncompressibleBackend::dsdTatPx(void){
if (!_dsdTatPx) _dsdTatPx = calc_dsdTatPx(_T,_p,_fractions[0]);
return _dsdTatPx;
}
double IncompressibleBackend::dhdTatPx(void){
if (!_dhdTatPx) _dhdTatPx = calc_dhdTatPx(_T,_p,_fractions[0]);
return _dhdTatPx;
}
double IncompressibleBackend::dsdTatPxdT(void){
if (!_dsdTatPxdT) _dsdTatPxdT = calc_dsdTatPxdT(_T,_p,_fractions[0]);
return _dsdTatPxdT;
}
double IncompressibleBackend::dhdTatPxdT(void){
if (!_dhdTatPxdT) _dhdTatPxdT = calc_dhdTatPxdT(_T,_p,_fractions[0]);
return _dhdTatPxdT;
}
double IncompressibleBackend::dsdpatTx(void){
if (!_dsdpatTx) _dsdpatTx = calc_dsdpatTx(rhomass(),drhodTatPx());
return _dsdpatTx;
}
double IncompressibleBackend::dhdpatTx(void){
if (!_dhdpatTx) _dhdpatTx = calc_dhdpatTx(_T,rhomass(),drhodTatPx());
return _dhdpatTx;
}
/// Return the temperature in K
double IncompressibleBackend::T_ref(void){
if (!_T_ref) throw ValueError("Reference temperature is not set");
@@ -442,6 +479,23 @@ CoolPropDbl IncompressibleBackend::raw_calc_smass(double T, double p, double x){
return calc_dsdTatPxdT(T,p,x) + p * calc_dsdpatTx( fluid->rho(T, p, x),calc_drhodTatPx(T,p,x));
};
/// Calculate the first partial derivative for the desired derivative
CoolPropDbl IncompressibleBackend::calc_first_partial_deriv(parameters Of, parameters Wrt, parameters Constant){
// TODO: Can this be accelerated?
if ( (Of==iDmass) && (Wrt==iT) && (Constant==iP) ) return drhodTatPx();
if ( (Of==iSmass) && (Wrt==iT) && (Constant==iP) ) return dsdTatPx();
if ( (Of==iHmass) && (Wrt==iT) && (Constant==iP) ) return dhdTatPx();
//if ( (Of==iHmass) && (Wrt==iP) && (Constant==iT) ) return dsdTatPxdT();
//if ( (Of==iHmass) && (Wrt==iP) && (Constant==iT) ) return dhdTatPxdT();
if ( (Of==iSmass) && (Wrt==iP) && (Constant==iT) ) return dsdpatTx();
if ( (Of==iHmass) && (Wrt==iP) && (Constant==iT) ) return dhdpatTx();
if ( (Of==iDmass) && (Wrt==iHmass) && (Constant==iP) ) return drhodTatPx()/dhdTatPx();
if ( (Of==iDmass) && (Wrt==iP) ) return 0.0; // incompressible!
throw ValueError("Incompressible fluids only support a limited subset of partial derivatives.");
}
/* Other useful derivatives
*/
/// Partial derivative of entropy with respect to pressure at constant temperature and composition

View File

@@ -26,6 +26,7 @@ protected:
/// Additional cached elements used for the partial derivatives
CachedElement _cmass, _hmass, _rhomass, _smass, _umass;
CachedElement _drhodTatPx, _dsdTatPx, _dhdTatPx, _dsdTatPxdT, _dhdTatPxdT, _dsdpatTx, _dhdpatTx;
IncompressibleFluid *fluid;
@@ -119,6 +120,14 @@ public:
/// Return the mass constant pressure specific heat in J/kg/K
double cmass(void);
double drhodTatPx(void);
double dsdTatPx(void);
double dhdTatPx(void);
double dsdTatPxdT(void);
double dhdTatPxdT(void);
double dsdpatTx(void);
double dhdpatTx(void);
/// Return the temperature in K
double T_ref(void);
/// Return the pressure in Pa
@@ -191,6 +200,9 @@ public:
protected:
/// Calculate the first partial derivative for the desired derivative
CoolPropDbl calc_first_partial_deriv(parameters Of, parameters Wrt, parameters Constant);
/* Below are direct calculations of the derivatives. Nothing
* special is going on, we simply use the polynomial class to
* derive the different functions with respect to temperature.

View File

@@ -857,10 +857,10 @@ void REFPROPMixtureBackend::calc_phase_envelope(const std::string &type)
PhaseEnvelope.Q.push_back(static_cast<double>(y > rho_molL));
isp = nc + 4;
SPLNVALdll(&isp, &iderv, &rho_molL, &y, &ierr, herr, errormessagelength);
PhaseEnvelope.hmolar_vap.push_back(y*1000);
PhaseEnvelope.hmolar_vap.push_back(y);
isp = nc + 5;
SPLNVALdll(&isp, &iderv, &rho_molL, &y, &ierr, herr, errormessagelength);
PhaseEnvelope.smolar_vap.push_back(y*1000);
PhaseEnvelope.smolar_vap.push_back(y);
}
}
CoolPropDbl REFPROPMixtureBackend::calc_cpmolar_idealgas(void)

View File

@@ -39,11 +39,20 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
std::size_t iL, iV, iclosest = 0;
CoolPropDbl hL = 0, hV = 0;
SimpleState closest_state;
if (
(is_mixture && PhaseEnvelopeRoutines::is_inside(phase_envelope, iP, _p, iHmolar, _hmolar, iclosest, closest_state))
||
(!is_mixture && pure_saturation.is_inside(iP, _p, iHmolar, _hmolar, iL, iV, hL, hV))
)
bool is_two_phase = false;
// Phase is imposed, use it
if (imposed_phase_index != iphase_not_imposed){
is_two_phase = (imposed_phase_index == iphase_twophase);
}
else{
if (is_mixture){
is_two_phase = PhaseEnvelopeRoutines::is_inside(phase_envelope, iP, _p, iHmolar, _hmolar, iclosest, closest_state);
}
else{
is_two_phase = pure_saturation.is_inside(iP, _p, iHmolar, _hmolar, iL, iV, hL, hV);
}
}
if ( is_two_phase )
{
using_single_phase_table = false;
_Q = (static_cast<double>(_hmolar)-hL)/(hV-hL);
@@ -52,6 +61,7 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
}
else{
cached_saturation_iL = iL; cached_saturation_iV = iV;
_phase = iphase_twophase;
}
}
else{
@@ -68,6 +78,8 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
if (!cell.valid()){throw ValueError(format("Cell is invalid and has no good neighbors for hmolar = %g, p= %g",val1,val2));}
}
}
// Recalculate the phase
recalculate_singlephase_phase();
}
}
break;
@@ -92,11 +104,20 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
CoolPropDbl zL = 0, zV = 0;
std::size_t iclosest = 0;
SimpleState closest_state;
if (
(is_mixture && PhaseEnvelopeRoutines::is_inside(phase_envelope, iP, _p, otherkey, otherval, iclosest, closest_state))
||
(!is_mixture && pure_saturation.is_inside(iP, _p, otherkey, otherval, iL, iV, zL, zV))
){
bool is_two_phase = false;
// Phase is imposed, use it
if (imposed_phase_index != iphase_not_imposed){
is_two_phase = (imposed_phase_index == iphase_twophase);
}
else{
if (is_mixture){
is_two_phase = PhaseEnvelopeRoutines::is_inside(phase_envelope, iP, _p, otherkey, otherval, iclosest, closest_state);
}
else{
is_two_phase = pure_saturation.is_inside(iP, _p, otherkey, otherval, iL, iV, zL, zV);
}
}
if ( is_two_phase ){
using_single_phase_table = false;
if (otherkey == iDmolar){
_Q = (1/otherval - 1/zL)/(1/zV - 1/zL);
@@ -110,6 +131,7 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
else{
cached_saturation_iL = iL; cached_saturation_iV = iV;
}
_phase = iphase_twophase;
}
else{
// Find and cache the indices i, j
@@ -127,6 +149,8 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
}
// Now find hmolar given P, X for X in Hmolar, Smolar, Umolar
invert_single_phase_x(single_phase_logph, dataset->coeffs_ph, otherkey, otherval, _p, cached_single_phase_i, cached_single_phase_j);
// Recalculate the phase
recalculate_singlephase_phase();
}
break;
}
@@ -142,7 +166,7 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
// Call again, but this time with molar units; S: [J/kg/K] * [kg/mol] -> [J/mol/K]
update(PSmolar_INPUTS, val1, val2*AS->molar_mass()); return;
}
case PT_INPUTS:{
case PT_INPUTS:{
_p = val1; _T = val2;
if (!single_phase_logpT.native_inputs_are_in_range(_T, _p)){
// Use the AbstractState instance
@@ -155,11 +179,20 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
std::size_t iL = 0, iV = 0, iclosest = 0;
CoolPropDbl TL = 0, TV = 0;
SimpleState closest_state;
if (
(is_mixture && PhaseEnvelopeRoutines::is_inside(phase_envelope, iP, _p, iT, _T, iclosest, closest_state))
||
(!is_mixture && pure_saturation.is_inside(iP, _p, iT, _T, iL, iV, TL, TV))
)
bool is_two_phase = false;
// Phase is imposed, use it
if (imposed_phase_index != iphase_not_imposed){
is_two_phase = (imposed_phase_index == iphase_twophase);
}
else{
if (is_mixture){
is_two_phase = PhaseEnvelopeRoutines::is_inside(phase_envelope, iP, _p, iT, _T, iclosest, closest_state);
}
else{
is_two_phase = pure_saturation.is_inside(iP, _p, iT, _T, iL, iV, TL, TV);
}
}
if ( is_two_phase )
{
using_single_phase_table = false;
throw ValueError(format("P,T with TTSE cannot be two-phase for now"));
@@ -199,6 +232,8 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
}
}
}
// Recalculate the phase
recalculate_singlephase_phase();
}
}
break;
@@ -225,11 +260,20 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
CoolPropDbl zL = 0, zV = 0;
std::size_t iclosest = 0;
SimpleState closest_state;
if (
(is_mixture && PhaseEnvelopeRoutines::is_inside(phase_envelope, iT, _T, otherkey, otherval, iclosest, closest_state))
||
(!is_mixture && pure_saturation.is_inside(iT, _T, otherkey, otherval, iL, iV, zL, zV))
)
bool is_two_phase = false;
// Phase is imposed, use it
if (imposed_phase_index != iphase_not_imposed){
is_two_phase = (imposed_phase_index == iphase_twophase);
}
else{
if (is_mixture){
is_two_phase = PhaseEnvelopeRoutines::is_inside(phase_envelope, iT, _T, otherkey, otherval, iclosest, closest_state);
}
else{
is_two_phase = pure_saturation.is_inside(iT, _T, otherkey, otherval, iL, iV, zL, zV);
}
}
if ( is_two_phase )
{
using_single_phase_table = false;
if (otherkey == iDmolar){
@@ -262,6 +306,8 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
}
// Now find the y variable (Dmolar or Smolar in this case)
invert_single_phase_y(single_phase_logpT, dataset->coeffs_pT, otherkey, otherval, _T, cached_single_phase_i, cached_single_phase_j);
// Recalculate the phase
recalculate_singlephase_phase();
}
break;
}
@@ -286,7 +332,7 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
}
}
_T = _Q*TV + (1-_Q)*TL;
cached_saturation_iL = iL; cached_saturation_iV = iV;
cached_saturation_iL = iL; cached_saturation_iV = iV; _phase = iphase_twophase;
}
break;
}
@@ -312,7 +358,7 @@ void CoolProp::BicubicBackend::update(CoolProp::input_pairs input_pair, double v
pure_saturation.is_inside(iT, _T, iQ, _Q, iL, iV, pL, pV);
}
_p = _Q*pV + (1-_Q)*pL;
cached_saturation_iL = iL; cached_saturation_iV = iV;
cached_saturation_iL = iL; cached_saturation_iV = iV; _phase = iphase_twophase;
}
break;
}

View File

@@ -63,11 +63,10 @@ A^{-1} = \left[ \begin{array}{*{16}c} 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0
typedef std::vector<std::vector<double> > mat;
class BicubicBackend : public TabularBackend
{
protected:
public:
/// Instantiator; base class loads or makes tables
BicubicBackend(shared_ptr<CoolProp::AbstractState> AS) : TabularBackend (AS){
BicubicBackend(shared_ptr<CoolProp::AbstractState> AS) : TabularBackend(AS){
imposed_phase_index = iphase_not_imposed;
// If a pure fluid, don't need to set fractions, go ahead and build
if (this->AS->get_mole_fractions().size() == 1){
check_tables();

View File

@@ -39,20 +39,37 @@ void CoolProp::TTSEBackend::update(CoolProp::input_pairs input_pair, double val1
using_single_phase_table = true; // Use the table !
std::size_t iL, iV;
CoolPropDbl hL = 0, hV = 0;
if (pure_saturation.is_inside(iP, _p, iHmolar, _hmolar, iL, iV, hL, hV)){
std::size_t iclosest = 0;
SimpleState closest_state;
bool is_two_phase = false;
// Phase is imposed, use it
if (imposed_phase_index != iphase_not_imposed){
is_two_phase = (imposed_phase_index == iphase_twophase);
}
else{
if (is_mixture){
is_two_phase = PhaseEnvelopeRoutines::is_inside(phase_envelope, iP, _p, iHmolar, _hmolar, iclosest, closest_state);
}
else{
is_two_phase = pure_saturation.is_inside(iP, _p, iHmolar, _hmolar, iL, iV, hL, hV);
}
}
if ( is_two_phase ){
using_single_phase_table = false;
_Q = (static_cast<double>(_hmolar)-hL)/(hV-hL);
if(!is_in_closed_range(0.0,1.0,static_cast<double>(_Q))){
throw ValueError("vapor quality is not in (0,1)");
}
else{
cached_saturation_iL = iL; cached_saturation_iV = iV;
cached_saturation_iL = iL; cached_saturation_iV = iV; _phase = iphase_twophase;
}
}
else{
// Find and cache the indices i, j
selected_table = SELECTED_PH_TABLE;
single_phase_logph.find_native_nearest_good_neighbor(_hmolar, _p, cached_single_phase_i, cached_single_phase_j);
// Recalculate the phase
recalculate_singlephase_phase();
}
}
break;
@@ -75,7 +92,23 @@ void CoolProp::TTSEBackend::update(CoolProp::input_pairs input_pair, double val1
using_single_phase_table = true; // Use the table (or first guess is that it is single-phase)!
std::size_t iL, iV;
CoolPropDbl zL = 0, zV = 0;
if (pure_saturation.is_inside(iP, _p, otherkey, otherval, iL, iV, zL, zV)){
std::size_t iclosest = 0;
SimpleState closest_state;
bool is_two_phase = false;
// Phase is imposed, use it
if (imposed_phase_index != iphase_not_imposed)
{
is_two_phase = (imposed_phase_index == iphase_twophase);
}
else{
if (is_mixture){
is_two_phase = PhaseEnvelopeRoutines::is_inside(phase_envelope, iP, _p, otherkey, otherval, iclosest, closest_state);
}
else{
is_two_phase = pure_saturation.is_inside(iP, _p, otherkey, otherval, iL, iV, zL, zV);
}
}
if (is_two_phase){
using_single_phase_table = false;
if (otherkey == iDmolar){
_Q = (1/otherval - 1/zL)/(1/zV - 1/zL);
@@ -87,7 +120,7 @@ void CoolProp::TTSEBackend::update(CoolProp::input_pairs input_pair, double val1
throw ValueError("vapor quality is not in (0,1)");
}
else{
cached_saturation_iL = iL; cached_saturation_iV = iV;
cached_saturation_iL = iL; cached_saturation_iV = iV; _phase = iphase_twophase;
}
}
else{
@@ -96,6 +129,8 @@ void CoolProp::TTSEBackend::update(CoolProp::input_pairs input_pair, double val1
single_phase_logph.find_nearest_neighbor(iP, _p, otherkey, otherval, cached_single_phase_i, cached_single_phase_j);
// Now find hmolar
invert_single_phase_x(single_phase_logph, otherkey, otherval, _p, cached_single_phase_i, cached_single_phase_j);
// Recalculate the phase
recalculate_singlephase_phase();
}
break;
}
@@ -131,6 +166,8 @@ void CoolProp::TTSEBackend::update(CoolProp::input_pairs input_pair, double val1
// Find and cache the indices i, j
selected_table = SELECTED_PT_TABLE;
single_phase_logpT.find_native_nearest_neighbor(_T, _p, cached_single_phase_i, cached_single_phase_j);
// Recalculate the phase
recalculate_singlephase_phase();
}
}
break;
@@ -155,7 +192,22 @@ void CoolProp::TTSEBackend::update(CoolProp::input_pairs input_pair, double val1
using_single_phase_table = true; // Use the table (or first guess is that it is single-phase)!
std::size_t iL, iV;
CoolPropDbl zL = 0, zV = 0;
if (pure_saturation.is_inside(iT, _T, otherkey, otherval, iL, iV, zL, zV)){
std::size_t iclosest = 0;
SimpleState closest_state;
bool is_two_phase = false;
// Phase is imposed, use it
if (imposed_phase_index != iphase_not_imposed){
is_two_phase = (imposed_phase_index == iphase_twophase);
}
else{
if (is_mixture){
is_two_phase = PhaseEnvelopeRoutines::is_inside(phase_envelope, iT, _T, otherkey, otherval, iclosest, closest_state);
}
else{
is_two_phase = pure_saturation.is_inside(iT, _T, otherkey, otherval, iL, iV, zL, zV);
}
}
if ( is_two_phase ){
using_single_phase_table = false;
if (otherkey == iDmolar){
_Q = (1/otherval - 1/zL)/(1/zV - 1/zL);
@@ -167,7 +219,7 @@ void CoolProp::TTSEBackend::update(CoolProp::input_pairs input_pair, double val1
throw ValueError("vapor quality is not in (0,1)");
}
else{
cached_saturation_iL = iL; cached_saturation_iV = iV;
cached_saturation_iL = iL; cached_saturation_iV = iV; _phase = iphase_twophase;
}
}
else{
@@ -176,6 +228,8 @@ void CoolProp::TTSEBackend::update(CoolProp::input_pairs input_pair, double val1
single_phase_logpT.find_nearest_neighbor(iT, _T, otherkey, otherval, cached_single_phase_i, cached_single_phase_j);
// Now find hmolar
invert_single_phase_y(single_phase_logpT, otherkey, otherval, _T, cached_single_phase_i, cached_single_phase_j);
// Recalculate the phase
recalculate_singlephase_phase();
}
break;
}
@@ -200,7 +254,7 @@ void CoolProp::TTSEBackend::update(CoolProp::input_pairs input_pair, double val1
}
}
_T = _Q*TV + (1-_Q)*TL;
cached_saturation_iL = iL; cached_saturation_iV = iV;
cached_saturation_iL = iL; cached_saturation_iV = iV; _phase = iphase_twophase;
}
break;
}
@@ -226,7 +280,7 @@ void CoolProp::TTSEBackend::update(CoolProp::input_pairs input_pair, double val1
pure_saturation.is_inside(iT, _T, iQ, _Q, iL, iV, pL, pV);
}
_p = _Q*pV + (1-_Q)*pL;
cached_saturation_iL = iL; cached_saturation_iV = iV;
cached_saturation_iL = iL; cached_saturation_iV = iV; _phase = iphase_twophase;
}
break;
}

View File

@@ -12,6 +12,7 @@ class TTSEBackend : public TabularBackend
std::string backend_name(void){return "TTSEBackend";}
/// Instantiator; base class loads or makes tables
TTSEBackend(shared_ptr<CoolProp::AbstractState> AS) : TabularBackend (AS) {
imposed_phase_index = iphase_not_imposed;
// If a pure fluid, don't need to set fractions, go ahead and build
if (this->AS->get_mole_fractions().size() == 1){
check_tables();

View File

@@ -706,6 +706,7 @@ public:
class TabularBackend : public AbstractState
{
protected:
phases imposed_phase_index;
bool tables_loaded, using_single_phase_table, is_mixture;
enum selected_table_options{SELECTED_NO_TABLE=0, SELECTED_PH_TABLE, SELECTED_PT_TABLE};
selected_table_options selected_table;
@@ -762,6 +763,42 @@ class TabularBackend : public AbstractState
}
TabularDataSet * dataset;
void recalculate_singlephase_phase()
{
if (p() > p_critical()){
if (T() > T_critical()){
_phase = iphase_supercritical;
}
else{
_phase = iphase_supercritical_liquid;
}
}
else{
if (T() > T_critical()){
_phase = iphase_supercritical_gas;
}
else{
// Liquid or vapor
if (rhomolar() > rhomolar_critical()){
_phase = iphase_liquid;
}
else{
_phase = iphase_gas;
}
}
}
}
/** \brief Specify the phase - this phase will always be used in calculations
*
* @param phase_index The index from CoolProp::phases
*/
void calc_specify_phase(phases phase_index){ imposed_phase_index = phase_index; };
/**\brief Unspecify the phase - the phase is no longer imposed, different solvers can do as they like
*/
void calc_unspecify_phase(){ imposed_phase_index = iphase_not_imposed; };
phases calc_phase(void){ return _phase; }
CoolPropDbl calc_T_critical(void){return this->AS->T_critical();};
CoolPropDbl calc_Ttriple(void){return this->AS->Ttriple();};
CoolPropDbl calc_p_triple(void){return this->AS->p_triple();};

View File

@@ -1,9 +1,15 @@
%module CoolProp
// Hack the AbstractState to return the mole fractions as a vector<double> which SWIG can wrap
%ignore CoolProp::AbstractState::mole_fractions_liquid();
%ignore CoolProp::AbstractState::mole_fractions_vapor();
%rename (mole_fractions_liquid) CoolProp::AbstractState::mole_fractions_liquid_double();
%rename (mole_fractions_vapor) CoolProp::AbstractState::mole_fractions_vapor_double();
%ignore CoolProp::AbstractState::set_mole_fractions(const std::vector<CoolPropDbl> &);
%ignore CoolProp::AbstractState::set_mass_fractions(const std::vector<CoolPropDbl> &);
%ignore CoolProp::AbstractState::set_volu_fractions(const std::vector<CoolPropDbl> &);
%ignore CoolProp::set_config_json(rapidjson::Document &);
%ignore CoolProp::get_config_as_json(rapidjson::Document &);
@@ -11,11 +17,9 @@
%include "std_vector.i" // This allows for the use of STL vectors natively(ish)
%include "exception.i" //
// Instantiate templates used by example
namespace std {
%template(LongDoubleVector) vector<long double>;
%template(DoubleVector) vector<double>;
}
// Instantiate templates used
%template(DoubleVector) std::vector<double>;
%template(StringVector) std::vector<std::string>;
%apply double { CoolPropDbl };