PlanckEinstein and PlanckEinstein2 classes replaced with PlanckEinsteinGeneralized. Docs to follow.

Aly-Lee not yet removed.  Next step.

Signed-off-by: Ian Bell <ian.h.bell@gmail.com>
This commit is contained in:
Ian Bell
2014-06-03 12:04:16 +02:00
parent 6db102b47f
commit bff7e2afeb
11 changed files with 396 additions and 286 deletions

View File

@@ -183,13 +183,16 @@
"c": [
0.6666666666666666
],
"d": [
1
],
"n": [
-0.197938904
],
"t": [
87.31279
-87.31279
],
"type": "IdealGasHelmholtzPlanckEinstein2"
"type": "IdealGasHelmholtzPlanckEinsteinGeneralized"
}
],
"alphar": [

View File

@@ -174,15 +174,18 @@
},
{
"c": [
1
],
"d": [
-1
],
"n": [
1.012767
],
"t": [
8.9057501
-8.9057501
],
"type": "IdealGasHelmholtzPlanckEinstein2"
"type": "IdealGasHelmholtzPlanckEinsteinGeneralized"
},
{
"a1": 0,

View File

@@ -193,13 +193,16 @@
"c": [
0.6666666666666666
],
"d": [
1
],
"n": [
0.944365
],
"t": [
74.91476960299131
-74.91476960299131
],
"type": "IdealGasHelmholtzPlanckEinstein2"
"type": "IdealGasHelmholtzPlanckEinsteinGeneralized"
}
],
"alphar": [

View File

@@ -168,7 +168,7 @@
"type": "IdealGasHelmholtzLogTau"
},
{
"T0": 309.213623675,
"T0": 298,
"Tc": 469.7,
"c": [
4,
@@ -180,7 +180,7 @@
"type": "IdealGasHelmholtzCP0AlyLee"
},
{
"T0": 309.213623675,
"T0": 298,
"Tc": 469.7,
"c": [
0,

View File

@@ -12,7 +12,28 @@
"level": 1,
"metadata": {},
"source": [
"Conversion of Aly-Lee form to Planck-Einstein"
"Derivation and Conversion of Planck-Einstein Terms"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Start with a generalized Helmholtz energy from \n",
"$$\\alpha^0 = a_k\\log\\left[c_k + d_k\\exp\\left(\\frac{b_k\\tau}{T_c}\\right)\\right]$$\n",
"\n",
"The ideal-gas specific heat can be obtained from \n",
"$$\\frac{c_p^0}{R} = -\\tau^2\\alpha^0_{\\tau\\tau}$$\n",
"\n",
"After some math shown below (and some manual grouping):\n",
"$$\\frac{c_p^0}{R} = - \\frac{a_{k} (b_{k}\\tau/T_c)^{2} c_k d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}}{ \\left(c_{k} + d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}\\right)^2} $$\n",
"\n",
"Two important identities:\n",
"$$ \\left(\\frac{1}{\\cosh x}\\right) = \\frac{2e^{-x}}{1+e^{-2x}}$$\n",
"$$ \\left(\\frac{1}{\\sinh x}\\right) = \\frac{2e^{-x}}{1-e^{-2x}}$$\n",
"which yield\n",
"$$ \\left(\\frac{x}{\\cosh x}\\right)^2 = \\frac{4x^2e^{-2x}}{(1+e^{-2x})^2} = \\frac{(-2x)^2e^{-2x}}{(1+e^{-2x})^2}$$\n",
"$$ \\left(\\frac{x}{\\sinh x}\\right)^2 = \\frac{4x^2e^{-2x}}{(1-e^{-2x})^2} = \\frac{(-2x)^2e^{-2x}}{(1-e^{-2x})^2}$$"
]
},
{
@@ -35,58 +56,155 @@
]
}
],
"prompt_number": 2
"prompt_number": 4
},
{
"cell_type": "code",
"collapsed": false,
"input": [
"a_k,b_k,c_k,d_k,tau,T_c = symbols('a_k,b_k,c_k,d_k,tau,T_c',real=True)\n",
"alpha0 = a_k*log(c_k + d_k*exp(b_k*tau/T_c))\n",
"display(diff(alpha0, tau, 2)*(-tau**2).together().ratsimp())"
"print 'No derivatives'\n",
"display(alpha0)\n",
"print 'First partial w.r.t. tau'\n",
"display(diff(alpha0, tau, 1).simplify())\n",
"print 'Second partial w.r.t. tau'\n",
"display(diff(alpha0, tau, 2).simplify())\n",
"print latex(diff(alpha0, tau, 2).simplify()*(-tau**2))\n",
"print 'Third partial w.r.t. tau'\n",
"display(diff(alpha0, tau, 3).simplify())"
],
"language": "python",
"metadata": {},
"outputs": [
{
"output_type": "stream",
"stream": "stdout",
"text": [
"No derivatives\n"
]
},
{
"latex": [
"$$- \\frac{a_{k} b_{k}^{2} d_{k} \\tau^{2} e^{\\frac{b_{k} \\tau}{T_{c}}}}{T_{c}^{2} \\left(c_{k} + d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}\\right)} \\left(- \\frac{d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}}{c_{k} + d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}} + 1\\right)$$"
"$$a_{k} \\log{\\left (c_{k} + d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}} \\right )}$$"
],
"metadata": {},
"output_type": "display_data",
"png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAABMBAMAAACPE1EUAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEM3dMiKriXZE77tm\nmVQzv5s6AAAJ0klEQVRoBc1ba4gkVxU+VdPv6ReKjySu3faaB7JxJ+al/rGJBGRRpjEkQv5MIxJH\nMWw7KruIyzYxrBuzy7aJIEZlRpOwSVzdSX6EIK7TRpQlRLeNiQYUpkB8EcjM+txNdjOe+6p7btWt\nqumanuj9Ufec75zz3dt1q29V9TcDsNVWfGirDLx+QjTjzWVpvPSo7AnRRNHb8GzTho6NTYhmrHGr\nu57n+Zl3PXiNrbDct6EhLImGFnyHOluwK8OPQq4JUC0OrCx/tKIhMImGFkz3qZfedjsvQXWE9fm6\njSS7yW9xAo1B7TxquKkdFz4FFQ/LXSvFVN8Kh8AEGjP/gY7pR3kHTllOZ/73v5b5xf0eTF3yS4DL\nrAQLVjQMJtCYBZWm6Ud4+Q67koPtC/BxApV+/jbiUTN7nnoJdjRNoNB5LQDY3alO5l/hyPUwSxbt\nvkPhDIFML0dFLHg0TTD5i0HA6pcGtjP3KMwNdPoR9+3aMay1tuHGO9E0wbpGP4jY/epZG75GvgEv\nlO+xpSDWisCtcDRNMD1vnRLNuuxB5olNovijEQ3Be6W3QZqRwJzMhRCkgW/9QNtAWDYIHGFm43hF\nzRWsOyDsmid6cczrD1K+8uDvrqQxbefWtR22zC9bDE24NHk1GXlZTnKxThlwN1StAnvgrcox+0bX\n9A2vaK56DI1Rxp1jnTBmIJz8G/B1Dn6Ohoq9Ylv5VfgFXKIcsz8zMH3Dy88YbgyNkcedRi+MGQgj\nz7z75cc5eDUN3fTyS3Xtz2vTtFZM1/TckelDJE0gD113KYwhkv3k3deIgDt/1CttbKxzTF6ZBfxC\nnYOTG/hlcvYf7fDE66w8CP49KgBw2/5ZT0STacIs+X+GMUQe6EyvQ+E4Wo13gLwdIXZWYPtu/sSt\nYsJwfyfbYwwOfrbql0bMNFvxP6ZPvOLDsFKHr7L6RBpSpszqOWXR3jkBlWWA9yGGe7cYnGFNjmU8\nwCBv2R/Cjjazik2APwObRaDlox8DLh3Cc7jvDnBBE2kCrMx1XrWAUDnLdnKHff2fA0dMnWE9gUGm\nK6tyJw7/lJsLNwK8U4JG55p7CI2t1gGrGnWAZBpap+zHlEH70jI02lCcQexqKL7GF5VhA4FBfiiz\na11psO5y6BBPmjVGYm/fBweX5H4MJtPYGPZahgPcd9agnh99qAsX+BIM8OxoDKbakqrWI5wfOT0k\nnjTnlsOYQHA1q+sOHMgcgmQaG8nBoQVtePBEteP+6mtn4TzgJckWlWAwp2pKOPW+cqz92sgKI4iP\nrdPNP8GuvzXqyTQ2kjM4eKi5/fINeaixl+XPwCmxqASDL6uKQhN2dJRj7Y/1rTADr4KfdXtwsQ2Q\nTGMjWexaUGf+87fdA3vu9QB2HBriIwwuKsHIjePoEQzHtdXo+Ju/cvP8IPPIT7A8kcY2xFzThnLs\n3sxSWwTZogpLY8JPPK4M4lMKy3sy8RnhqCsmU1sKhyTyVGamJ0y2qKJpTCEJ/e5OfMJ0/4P5+IxQ\nNLsip74eCingcphvc1ssKjd9TCUl9SfFMJFpbid3d2TQDLxJuM6pvYKzFH3L8AvTLKoq3lDG1vu/\nKopjYupTlhdmlaL6FIuqSuF/PPUxFtWfsjK2cequ/dFRjbzV3tnOqeNjxDY255XJkQev9enz8IZr\nWcPnPvpKPgGbTVpO3dnJhnhPO90Q8tNbpj658xJmmtgFs9Bqvb/VEm+YcofBs76t7f//Wo+UsNJN\n3UoXvGDid5jycLNLskQSna520k0dKJ0iC059KnZzvEuVJfXZJs24TztJDwI6k1omnYyoqf/25A1D\nBsXeTYvflVWJHZOwck2Amx67Fu/Sha5f8AFx0/b9zRmKzhDW1NQVRWldWZb+0rYFtEJMwmJykvdN\nHt7lJ620fXMMQ9EZwtq3AwS1uKnfGUiOdpmEVfEw/lmeM+ef61UGjt18OruwJvgaM9G8mU08Vcpq\nJmFxOanJgVxX4hDzgqdSLL1P51qCCppdVla4D/0iGE5RCJOwmA7kCDr9odd6KmWcXtFFCGuCaq0b\nplTa3FwnHItBqA50lcqb7Spr3J7SWWttC6q0ud9YKyJBqgM9o7IaTWWN21M6a+1qOwwrbe6RcCgO\noTqQr5BNrceVxMUonTVvN0Gzrzy599U7dvekNpeNvVsFhSXCgyb7UYi3SiyJTGKdISsRPNr8Bwnl\ne1CagcqQ/76LMkz8Cat5pJQ+IzOY3Z54q16QRmInf7zneUE6W3HmIkFzAHMjKA9AaHOVGRKTptPx\nscW6bzLj6X0XF3o+UliXZvbfPmYxCF1AVkpWxAr0s34Y4IwHhToIbW5qOTxYte9jhrAEMMosOTpY\n9d/WuXjpFwUMQseUH9KSFTFDGhwA7K1DEUBoc7URoZImGcsQlvARoNx1kEC2jH9KDnYUZukJXeAm\nkqyIlZomIT9HUpur9Xls3+G2zlFjodikJveWwy/weK6v08DxL/G5HoHRtNLZZKV5sy7kLfYNKMsv\nfanNNYYsdovnNrmqxBPV1LWwVLwePsZD/o/u3DvHj3hwZ5TFezudRVaC64y6sIMXCG1lfiq5Nodf\n2SELHYecx1UlnienroUluN1znuehBh7/Ii41tISKg0bZ3B2tdOw3/KCsFCms8cHwQDcYdHP0gWvW\nQwQlJHyv92E5dSIsffr0EfHx8f2i2vV/+tTE/odAogg6WA3LSsVmhLDGaLAFTgnb1nWbHaJdWsZD\nUcD5Vmvn91qtEUeVsKQ2v/zOHrhtTBZNT/1gXWHY2+nAIitFCmuSzW0SWjRxW9eNXzBczuaqEg/I\ns06EJXpS3U5WlWvYv7GykJ0uSlayCmtyhMW2NGSH27pu/I9lakPc97iqxANq6lpsOo7vc6om84cj\nyvR3GCj7VxvG7HRRspJVWJMjBPZmph7jn/wd7vNwrYddpYuiLleVOCanToSlu8ARmyOPq0NWbZ0I\n3KhA7O10KWSlwB8X5Z/YuKKPgnRmiQ9W6rJu37N9oSpxTE6dCEuFF5/lEfOQOa/92Y627XQwvqxU\nGhJSZVaWnTa3p5sK0gqSnLoIaFgl+j09KeUZHw7WGXRjKRB3BDi5W/Ikml+SBmgFqThUGPYaJqAw\ndS36dwbCus6gG0eBqAZPBx+iptZXn7kIBSkCZizTXc4lDq5HHDQj6sZRIG5XkzSY3QGMOKCfQ4z4\nphy2k/hNXoC+PwHjVitH9sXTdRF40hrfFGh8NTdVMdmkZ9LTLaQvnUgl/kVD2nYibeGE6grNtETV\npbSVk6p7KC3R1DBt5aTqnk5LlLow7YChurT/EZF9OET1egNO0ktWxITe6EUEXkf4FutNK3ECP07M\n2P6EbDvNGORHjTTlm6/5L9dwOhTAi9ZrAAAAAElFTkSuQmCC\n",
"png": "iVBORw0KGgoAAAANSUhEUgAAALYAAAAoBAMAAAChsh0EAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAIquJdjLdEETvu2aZ\nVM0GsGrEAAADx0lEQVRIDa1WXYgbVRg9+ZlMNr+DSosvZqBFEKQGd6tvy0C3IC2WEUToixsVrb6l\nFAoFIQERLYXullZE1ocgfRFR8tAX6Q8Rf0AWYUARH0RTLPqgLLVr2UK3rN93586dO9OZTHbpB7lz\nvnPOd/Zm7iRZ4EHVafdBJd2fM3M/lcVUh1kOqT+m+0xb79LwkTQhzv92xb8pL3/0Rhc4G5cT+vKH\nCWQitZobYx8p3vssF21eJ9epqU/oen2AXznsWZH49+RcVp/MtkjHON/HD4zHgmjRjZlcxsZkXVPf\neR44dmYEcyDIvK1pibAxSqTTyJXdw0DK3lXLCrzTXMs37NB2I4TJ6LtkOoUtDv/ylPSlQingQgqf\nTNesS32ltELInLGvrSQBcrejfbQrfzoKieZWWDeZnollnbRDM6OqcEU5rWt6WgNcm92cDwPz44iI\nih3ta51oz93xkFrqhpjQyFg2HcUUY/uKZ5cGyqrATwrhmRAy8qq22VdU/V8FBYhnN0dRnTst+5+Y\nmnc0wrgH5PYcPBZQlG3OHeijvGd1YUBk0xHK7IIrrmIJsmmOpkWdWPhKXEuaDeZd4LRVuYnieSFS\n9ksWLmKpWxo6xPSGTO/yGmNpoC7IprkNf678Ch5mH3r0+vGKgLT8B3MFtQHwpmAo+zVg0X0c+Q4T\nLZF9HnlPGoiT2Tw39uce9cxv2c1f2nW7ICAt66htoOTA9L+TKjbfpWb7e5m96JGHHAgMhGQ2z7V9\n+olLh7skoPBnGw2XkahNzAzQc1HuiLZi1+9Q9qDlnuJULA5pIQcCw/za2ltra3yGPEcH0yF4i15B\nNaxcADfRa+N3dAujF23i5L5Hj8ztFw5xT8hBe/INhOS+tbl14fUXY/Vw0K2j5+GLutW4+h6/84qN\nj4Gl/i9S7zkEmrR5TxqoDbLDOXoOinJAv9xFw6m+WkDTP4zKAEdd/IGz774gXE3ecs3G2640UCuz\ntbkjMP1nUE9G7h7Mvfsfeg5HD3kkGE/dsnJfP91HY2uL9k831eZ19rIjDdzJbG2u+M1lFmJl0MmJ\nOmQsu5r2OY7/7FBfGQdkaJDZvhDSgVFd6/6jB1w0Ovz2ZdEvNqr0QmFZMprhTEDxNTqnK+HsdezV\n9m2+Tlv2yKn+NqIGFZJCs07PRWKdmDsgDpe/E3Za/HhNrE8mqhPFRWuiDGT+oKbPz6dLvkKfrJ3W\nStZgcZzlSNPry2mK4j9QaJuglHWU9Nu9zUhln2Kw4ij3tkDuXLadP0Y7qZPTPAS7rJ1E47NppnLu\nNK64R/yf8j+av+zHorG6yAAAAABJRU5ErkJggg==\n",
"text": [
" \u239b b_k\u22c5\u03c4 \u239e \n",
" \u239c \u2500\u2500\u2500\u2500\u2500 \u239f b_k\u22c5\u03c4 \n",
" \u239c T_c \u239f \u2500\u2500\u2500\u2500\u2500 \n",
" 2 2 \u239c d_k\u22c5\u212f \u239f T_c \n",
"-a_k\u22c5b_k \u22c5d_k\u22c5\u03c4 \u22c5\u239c- \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500 + 1\u239f\u22c5\u212f \n",
" \u239c b_k\u22c5\u03c4 \u239f \n",
" \u239c \u2500\u2500\u2500\u2500\u2500 \u239f \n",
" \u239c T_c \u239f \n",
" \u239d c_k + d_k\u22c5\u212f \u23a0 \n",
"\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" \u239b b_k\u22c5\u03c4\u239e \n",
" \u239c \u2500\u2500\u2500\u2500\u2500\u239f \n",
" 2 \u239c T_c \u239f \n",
" T_c \u22c5\u239dc_k + d_k\u22c5\u212f \u23a0 "
" \u239b b_k\u22c5\u03c4\u239e\n",
" \u239c \u2500\u2500\u2500\u2500\u2500\u239f\n",
" \u239c T_c \u239f\n",
"a_k\u22c5log\u239dc_k + d_k\u22c5\u212f \u23a0"
]
},
{
"output_type": "stream",
"stream": "stdout",
"text": [
"First partial w.r.t. tau\n"
]
},
{
"latex": [
"$$\\frac{a_{k} b_{k} d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}}{T_{c} \\left(c_{k} + d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}\\right)}$$"
],
"metadata": {},
"output_type": "display_data",
"png": "iVBORw0KGgoAAAANSUhEUgAAAJsAAABMBAMAAACIdNgoAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAIquJdjLdEETvu2aZ\nVM0GsGrEAAAFE0lEQVRYCe1YXWgcVRT+Zncns5vd2QwFfTSL5kHwoQuJiChxiSnUJsgKIgTFrBZ/\n8GkQKbVGdp+q9aFsrRY0KkH6pqarFEVa2xV/XoS66ov4ku2PFUmJKRFSbUI9587O7J3snZ+tPvTB\nA7nzne989+yd2bsz3wToO9Jv9j0ldEI9tNpvUWv3OyNUby5/I+qPv/W0HSqMV8w1bkGiDbRejaeP\nUOWtX2A2SXRXhDBeOY/bkWuRth1PH6FKj7UwsPckjIUIYfxy5uOX4oujlQcnojV9KKbyL/ehjpR+\nmd0ZqfEJ9h2R0tfekZJrUkh0BDwn1zflBKdGN8aLPiY6kTuk1336pl43Sj4mMvF1SFZ8+la2YNR8\nTGSSrEiSfFNKCCZKNMQI7bYdzziy/Mg0/5A4to3NdqAxNm0xM1DmMToOWIOrSB0mYfU3XMYDRULp\nQzhmO3C/pTGDKv39dJJRaBjzyC0Az5JoycYaBmqEXmzgDC2IoPYu5sS6DgJmIUm18MitY6AEg7/H\nMzDWULUJfW9jNy2IYGJ+8jNukLxYRF705Sw4MguolpGukOIS0pvYz9K3YVyBgEMFzp3IW5oLA4/V\nIpZgJ5sPF3AVtNR79QlepLlqCDhU7M7Uv5vqJgGo2sKiaeW/eGUdV/iaLf9ctY1NDLZ/EDBD7UoB\nU1V0vpR9MokhfrTcgQ+BDb5AF/BpoShgqo05SzUvgDNG7ty2EzO7WsDcRAP6Gx+R8Pn7HxqpOXB6\nqhEwM4TepdfLopxamNE7OgmGzFSWjuqVoigMlh5xN5cElXNCyLMYcVaXtxI7OjoJhsz8v3QDXQHp\nsfTv4eoNdGL/6VJiOu4TcT+0LgmNgpTIkG5C8cLvuOnpoYwDtpLuJdlxJ9qA/vsRepamCr0KYozY\n7yDsuNkim+kad1rmoScSzR4qgGDHLSxyUpzPsPKsZpWsqiM7brbIyItqoqASXVKRSo4dt7DI+0RZ\n95ssZ4p2WTk1iJQt8gWFKFdRkMGUbJE/V8gyJQUZTMkWebjGOu2v97f/fesfRWfOcNk5qkefRd4q\n4Wc62ZciMhXkGk719FaRP5cNrnw7YxVvaRqB4SayNcbAt85BHg3Ly3wGF8iev++r814xtcrwUeB0\nCynboZ/zqh4wSx5MVjzIIIcZvOAxptgUNWC7jXSH3fCqHpDabbHIJj7BXk+nu1dCeiH5U1RHJ8ue\nCl67XouMka4MxlUn0bpLcqibW/m28Mei7rZTWGQ8JbWDsxRk3VXSrhGfcBiJlvDHvnZk985sscgG\nXX3znman55pzTHR/beL8yXXSjcoj3dX1WmSk28CPKHbadU6Stp0bOq+ODDK5/woNtCtXVi6+t7LS\nJKiwyOO7gV9ZJqLTjrYdzZss0ShOlgwy5eyPRXRWF2CRcRaWo6P3BhG07ehVQq/TaPDlHGrQPzWE\nPxZlt12ARX7sOMlFdL7ZYzbtyAWjzBx/Qq5ALyPCHwuVe+2iLLImvtHk4rVzJWRaYipe58PoiZLj\njwXntouyyDq9fLgxZDlo0SU8f9zdxlwKscimtxfofl9Dk+V388Dh+WOkGw4jxhCLnKx3ddrXx23O\nlsRIwPPHXQ2jEIs8WPBLOavKS+kthzG8IbaG6iO2agLyWau3kK33cjGZcYWO3ueuN+ZVEx9UkXE4\nU3lemXKcuQrNQENBIttWsTG4U2rNUTUdxWqH1Io913e2e1rqdtoTaj6C/SCoflNQIYw3ak71H3/h\nm+3fsAlYAAAAAElFTkSuQmCC\n",
"text": [
" b_k\u22c5\u03c4 \n",
" \u2500\u2500\u2500\u2500\u2500 \n",
" T_c \n",
" a_k\u22c5b_k\u22c5d_k\u22c5\u212f \n",
"\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" \u239b b_k\u22c5\u03c4\u239e\n",
" \u239c \u2500\u2500\u2500\u2500\u2500\u239f\n",
" \u239c T_c \u239f\n",
"T_c\u22c5\u239dc_k + d_k\u22c5\u212f \u23a0"
]
},
{
"output_type": "stream",
"stream": "stdout",
"text": [
"Second partial w.r.t. tau\n"
]
},
{
"latex": [
"$$\\frac{a_{k} b_{k}^{2} c_{k} d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}}{T_{c}^{2} \\left(c_{k} + d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}\\right)^{2}}$$"
],
"metadata": {},
"output_type": "display_data",
"png": "iVBORw0KGgoAAAANSUhEUgAAAKgAAABQBAMAAAB/pCPzAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAIquJdjLdEETvu2aZ\nVM0GsGrEAAAFkElEQVRYCc1YXYgbVRQ+yWx2JptkMhQUfHGXug8FwUa6PoiLhnULtV00ggiL4qb+\noCjCUETaurIBsVIfZLVaxPqzlIIvdRul6EMrRmwVfKjRJy2UTVt/kNZ1RWGr7RLPuXfuzNyZzJ2Z\n6EPPw733fOecL/fMT3K/APRtxht9lyoK5xWxfkOZTr+VirrShZMsmvt1/+OKtHShYvN6yHYASkYj\nXaEq27ROQ6mFGZqtSksXM2EDFNtYY6aqGz+s2oMx1obBnccBdqch1Sy6ZErLf/icMh4ODlq5P8Ko\nhLwyIbkJnHwjcykmbZv5fExGOFxaDWMS8nlhi+Srnd37KT5YlbNeftvnd33mg1XLcxQcD2asSUDh\n/B1fnJeQGIfKC61AkiFfjiJMwzOBFJXLyl+AF+UcrS75JfgIdkqI2tHqALnfzrwvZ5kt2YfRgN/L\nzdyw2fnCMUen2vludwWzNk3WeO66sRl6K9H0sSmLLR5jo3rYaw2twMA+TJr7BZwH/9q22WGQ8Soc\nseGuCkb3WBmaQF/BL6nbWrSMNP0AFBcAnsCEJRv+5Hn7INtm0K4mnMKnrAGQeQdmaxQ1OgDfAeMn\nt6cVV+nJ1OkWnwKdkyKGO6LhGxsewRZsgOyByU8Ywe0I/MRW0UN+AeZqYNQx4yIYa6xVxHBHBL0F\nOr61e3BVHsFB2FmwxLLnPFeBJbC11n0jcBlo2w3cGTVHEO68tKLDeG4Cyv6GHzza7EkmwLk2LJYs\n89OXVuES4CWkVstY0iZIX4Ohzrdw4fs5O4+kVVEUN5vVwsMalE9i3o1wmLdaHMGbwqAf4eORClyp\nAQx0YNaKIxNxffSWdVtgemsbYHaiia8+tgqbjlU5tOPOe0cbudc/wOypbRhNZVtz8zVeQK2ylQcN\nLEznUrE5yYdy9QpfUqvMPGioer/mgKmmszDKyXirVOtCYFrZzanYgsn9thrkkfx+W5VIgs5/bjVI\neBX6vp/E/215FbbZY0vqQ12PggRQgkNdApZASoJDXaAigZvgUJeAJZQSe6gLVUhAb3UWPNQ5Jcek\nUoUz74/tcJzQoY7h+GORzGR1RmcAtNChjqGw1+Zz7EjqLNvBMxdXZ/jTgxY61DFUT6xjSZ2RkHLU\nmdnE+vChjpFmMS2ZkTpjQoqrswwdY5xDXZBgJmn3KMo2cCHlqLP1QSrPv+gtY1akzpiQctRZ9HYy\ncdJJ/iS/kDKrcgyPS19/xaBiPRhR+n4hVZCeWyq7GW5i1fkqm5IOfiEVlnLbYcYipuFaUr5wHh3B\nmWX+Prjxn/W/V96E4QYBnzE0alCrM94r1moVyNeh2GQncvS/jOLjuFKd3SNqs9hxCwoNgEcZ9KQI\neLPOrgvz1epsSSQ+gB23YcDG0zIru+KRiVWpKlag1d0lLQLqTCgiaABstMEAOM7T/+KTf/SRqtXZ\nHF5EYUwwGxWjhoB+maGuICPPJY1TZ+UqK6Yhwzp+6MxpmxxG6goyliRIY9WZT2AU2B19r9slhhxz\nXEEmke6KU2flFsunIUtfWcJytFNPkDFY7DRWnflI8TH1jLXvCTJ8jpeXf353ebmFKbHqzNc+Pqae\n6XT3XUHGcWen8eqMNJdj+Jj6jBQpBZkg47ggjVVnwx7TERtrtckqZ3gNJ0+QMUxc01h1NtPgHKAt\nds9VUY7n5jmwSJMryBgmSGPV2RJtz7Pigl7j3q0C9ASZIGURlWRxazlHvi243E/zBJnRFEGcVZLl\noC8Pl2VL+O776wkyEWKzSrLwV9xNNxvQ4s7QiAumXQT/bsycOGpzjvCvV2JurROViv8V9Gt5/7WX\nSe6W3RTeD9G5+Vp0TB15Kjpc6ETHlBFl4SFlaXTwOvex7JHzbH/960/34HKhzHZ3mWYxFH3vieaa\nNFxuru+P1H8Bv/Dea92wCuIAAAAASUVORK5CYII=\n",
"text": [
" b_k\u22c5\u03c4 \n",
" \u2500\u2500\u2500\u2500\u2500 \n",
" 2 T_c \n",
"a_k\u22c5b_k \u22c5c_k\u22c5d_k\u22c5\u212f \n",
"\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" 2\n",
" \u239b b_k\u22c5\u03c4\u239e \n",
" \u239c \u2500\u2500\u2500\u2500\u2500\u239f \n",
" 2 \u239c T_c \u239f \n",
"T_c \u22c5\u239dc_k + d_k\u22c5\u212f \u23a0 "
]
},
{
"output_type": "stream",
"stream": "stdout",
"text": [
"- \\frac{a_{k} b_{k}^{2} c_{k} d_{k} \\tau^{2} e^{\\frac{b_{k} \\tau}{T_{c}}}}{T_{c}^{2} \\left(c_{k} + d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}\\right)^{2}}\n",
"Third partial w.r.t. tau\n"
]
},
{
"latex": [
"$$\\frac{a_{k} b_{k}^{3} c_{k} d_{k} \\left(c_{k} - d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}\\right) e^{\\frac{b_{k} \\tau}{T_{c}}}}{T_{c}^{3} \\left(c_{k}^{3} + 3 c_{k}^{2} d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}} + 3 c_{k} d_{k}^{2} e^{\\frac{2 b_{k}}{T_{c}} \\tau} + d_{k}^{3} e^{\\frac{3 b_{k}}{T_{c}} \\tau}\\right)}$$"
],
"metadata": {},
"output_type": "display_data",
"png": "iVBORw0KGgoAAAANSUhEUgAAAZYAAABUBAMAAABTgQyuAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAIquJdjLdEETvu2aZ\nVM0GsGrEAAALOklEQVR4AeVafYxcVRU/s7Ozb2bnU4gaEmNH3KQaA4x0xbQ0ddKPBGhTxtRoGgy7\nVUFFsSMSQrE6kxirVVMWK+WjqBvSRIzYjqaRkCodhGJSSBmVP5SE7AC1xrSuWyUu2tb13HPvufe+\nN+/Nx9vZbUlvsvPOPb/f79x75915M/veD+A8tvj9wYO3w4JV5xGZaDN2O6yNrN/QoS4LRpo28Q67\nA27MBS1mZ7DZ5Wjpk0eIGfvbns8CxF0yxm5+6NZyl+UWgraz28FTtXfDAK4gHa+KeawWL9wYa3yb\nM+fj6LT5RLvnk8m9DOk65qK0+kzNgjX2YSu56OFAvdshM/A+SDWQnSFFZNYSaqxpJRciHL26zTYa\na4O5JxMfbcDQ3b8GuEfmL7dgxpxJK7kAYaqWajPCqZ5GTPziq5rvfRdsTJP6HKQKyZnAkpHTgZAf\nsMv6xGeKboaNuZF+9lLNwGqp8UDID1if+ZpOJyd0SIGNuZF+9raVAqslioGQH/B08jqdjrypQwps\nzI30p3fPHqxz53FPsfhP65xZErxMpHznEebhcc5qlN5tQAuaM9k+R6+Jei2f/WyDhznMgf/xnCud\nfP0jz7xuMleaEMCD2VC/YppMzD0lgHvLXP85DnyPcfs7BN8T2Ax3GuKNJmzBbKhPsZjMtlLkjKfc\nNbp/m478gui4K5uGX8LdJjOVMzF4MAvpVxgdBxjKJR/w1DNfKmc9iLubqbv7MGL3x/ROpawLs3nz\niiPvXYe/ZEXLjGxoxEavFm9gdHQdpQBRs+feoJzGJEO9XjLKs3VGN8hzcItNqNREzx+zefOKd+aG\nZ2BQXGYqfwX1ZRi5Habh+gLmEJ2VKE6Edp/B7GHj98GBspTsyEWEEpwZ/Jl8bV2xskUR+GOKMu+D\ns5euXJ/DQlNl+JeslyjAFTBUxfkg2gQQKID8JGmMcvyyrQbHcIdWkfVD2F4SafFPyx+AloW9rAgC\nMET60lKzMFQER1yDjoGj1vJiGbsVfBFoQaIA8gqnMVzpiVuxfamEZMx+WkoG9q59gma2ChN/oUi8\nZOv4EoAJvB8tMQmVEsTHsdYpiJ+T2+RhUXkH/gm0KlFcC+0xjQmObg+Dg9/sQpLN6yTAq5BTPVpL\nAGYJ5hVWCjAF5Wj9Y3k4Q6ehiu83nZ6VsdVgoWqPGcweFrPpGQeEhDYTY588WFMhpQMwZs/3WGnA\n/nQu85tvzcKbgLte7CxcFX4FnPxTpWyhuES6jmnMHtk5B8PN35MEP09QtDEZZ8WiArBWdrhMppj8\nVBSyR1D9AXhcbhP4LThH4GwJL9MGRQKdLo25xjsOv8oXSIK3N7bnXBh1ljTwEIC1ssNlnJEPXXId\nbL4Bh9q+uga0TSB67UGI/eDneCoMiuW/L4bQmGu8O9Z8dKRKEtiwHsu0tLGqSPljLeR5JW6ITZRk\nAbGzKBqc3BxTNRndz2NYGKfEMSBNlClZ1aYvULwvNl6QpcXOojZc/HhUhYwuV32wME6JY0CaKFpr\nCxYkfhVG5BrkNhFjZHID69RYjOr31sLs6QSkifKoTVycuN02AfmbKtxE6LspnDSsqt02geF82LIQ\n6+22R+hxbGG7bQLeGxC2sEMcbXYgLDaM34hhW6IWVrlQuo2hC/85tHKhhAm+YPc8wBd6Viy0INkM\nOUJoYcjxupHt64bkw7ks55M8z6m7wm0y5/bzPG+/4SNb/LIdc8MX3FVMTPntHeftR7Duk/nBF1nO\nvgf9Fo8vsjN3MS43vINl04pcL+KbVpd6oYc6FxOhVCiaxDtvPYgfSze7ph8KN6fwDpamM9uLuBwt\nuOlOPmjGeE8nVJuHgyU9zuJuRo6tAA99V5CsazuLp8A8HCxfz7HYU9O/G2966IN5f2L3dhaPPryD\nxVmKN0GUNcZT1K9bhNNe+kk/Ht5mr/vnO2bDO1i2v7OgxR2HAWciNuGlLyn76sb8075cVzK8g2Xr\n3DiLXRUDOi/hNZmdNIoykPflnlLZ6JqlvrhIrnw8YME9OFiiR3/nqd9W3Dodmx4Tj4Z0Y662szwP\nRzXoCaK5gaYnpbo9OFg+CFd6arQVt07HRT9u12KutrM8A4Gug6Fc0J2sHhwsW2AsZ08AgMSRrXNz\nVXeeeq3TcdGftCXMTRR19piOPEGi6nWuMKEHB8uDsKTKMnkkceoTSwP2r3c6LjrVivzn0av+e/k/\nCuJJpGjGztJycY4sW9GQo0Ka9qflgsG8/R+DovHh0mX4aMNupOQ7z1lLCWWYtIkyFnQ9HV+6eAAF\n0QIkxiFVY67eWM7GoqfouyB9RqWGJGZcMB6quxupw1TNnRLKz7hTqhcXs/K2bKN1OpKj6PQBHsBT\nUYdklbnPmTL82XdyMndjA/6t0JXyeG/ZsDF6atnZVa6ZKOXgOciOu5jCPxOt65Rt40nVOP1lDqTd\nhqfjdv0o+uAMkm8CONyAQZyU5N5mCow1ZJwuyuPOqvM/GSXr8mhcMNSvxyacokTkq1LGH4RK084D\noBJ9ltxsi89QibMvcUB04Om4XT+Knha31asAV5UhjpHksp0Ft+eSmizGa8FnrWqPfQO+SRB/F0ke\nNJJ5p6piOhild4/ho/VCXM/atvEM6QLWWk5Z0xGPfI0jSNHZOLUH1Zor7SzL1pa2wJNlWdbMSDhz\n0OkS+/srP8Gn38YF85W1TxN1oEgH/WKUYtVYVCKkvPmVl7E+K0e0xgS8FqLL6QTSpXUFInQq1NRl\n7h2NTHPTmveosnpGm76IK0AXTGJubsZ2wcS3wKXENftDSlkZu6YOIIoO7hYA+Wd+PId+Pa28RQpc\nr7wWotN02tDlOUjSQwY1dWln2Q0DDVOWZ4Tn5H7x9P0KgiwXzGUNfDouWgX//mh9DoxyeQmoqDDN\nKKUQsNJt8REINrWW7uj0KB4G7N8ytO/Q4GI3MyN4qPwibgxqlgvm/QfXy+wu3Mt5fhiLLKNMPEBF\nyVIjlVSElfGmbfEhSK+lO7r8nCfGlVYcyM6SmMRIbgeITk+f+NH0dB0zbwM40CCny/UFlwvmnwiK\nFj1RgExJxtjTykgOhs9RUfzHnpVUn5Vui4+osGp6+vPT0+JjhsJKVU4nmA5yLfj1YhrtMTS4YBPb\ngRq/u3NlOPAEnUthm7JcMPL8SnImF5GBeFXK7CwMn6GiZKmRSqpvK43FRxVQe6w7uqx0uGEGl3YW\nYUtp0HYghNeCp2RrDq/KafLDVIxHZjeeRa4Re349h3ot+M4OnaaiZKmRSqpvK43FRxXgtdBAnejy\n2+JAGbXRtUVZQawvlUenGm0HyvFaDkH6DemC2YF5ywWzERx5TZYlzKtSJmswVaSiZKmRSqofqBQ1\n1Fq6okfEBSy6f+61Im6o2AR2sJGdZdmhIv7EQIcVNV5L/IWjVXK6kB/GcsEMPntIMr2vrPzeC4/h\n1wsWJUuNVFL9QKWoxNcxstt0oMfQmKZaatIpyXA/p2g7UIdnxAiQbUr22AVjMDtqUVp0U99W2LFa\ni0x1oMsf7kRNNLjIcg5oO1AnXuMcH8+qhQOwC4YR97FFadFNfbfG9L5rQnQCym8wO2XH0Qndy+Y4\n5H8q5HbgrPto/DDoL1QeGTcjqGfRabsF8VrzHejDeS3JVKEuOxU+CdZ20DQVtPXDeMlB/Tb1/SQd\n6OIqqVrk2YNlGeoFWtuBWXxs64dhUqdjm/p+0g50770DKpGcUJWs7eCt3dYP4yUH9dvU95N0oK/y\n08zDzuJXbrFye30H2uibvcCTad5N7nmGt7O46yxqb8h89O1xk0279xaJnwqY576A/AWcjtwXMLmQ\ndpaAaouSvqsRMExIO0tAtUVJ/yxwlHB2lsByCw+oG1n/Bx1oNENKOewqAAAAAElFTkSuQmCC\n",
"text": [
" \u239b b_k\u22c5\u03c4\u239e b_k\u22c5\u03c4 \n",
" \u239c \u2500\u2500\u2500\u2500\u2500\u239f \u2500\u2500\u2500\u2500\u2500 \n",
" 3 \u239c T_c \u239f T_c \n",
" a_k\u22c5b_k \u22c5c_k\u22c5d_k\u22c5\u239dc_k - d_k\u22c5\u212f \u23a0\u22c5\u212f \n",
"\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\u2500\n",
" \u239b b_k\u22c5\u03c4 2\u22c5b_k\u22c5\u03c4 3\u22c5b_k\u22c5\u03c4\u239e\n",
" \u239c \u2500\u2500\u2500\u2500\u2500 \u2500\u2500\u2500\u2500\u2500\u2500\u2500 \u2500\u2500\u2500\u2500\u2500\u2500\u2500\u239f\n",
" 3 \u239c 3 2 T_c 2 T_c 3 T_c \u239f\n",
"T_c \u22c5\u239dc_k + 3\u22c5c_k \u22c5d_k\u22c5\u212f + 3\u22c5c_k\u22c5d_k \u22c5\u212f + d_k \u22c5\u212f \u23a0"
]
}
],
"prompt_number": 4
"prompt_number": 16
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"And then with some manual simplification:\n",
"$$\\frac{c_p^0}{R} = - \\frac{a_{k} (b_{k}\\tau/T_c)^{2} c_k d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}}{ \\left(c_{k} + d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}\\right)^2} $$\n",
"Aly-Lee starts as\n",
"$$\\frac{c_p^0}{R_u} = A + B\\left(\\frac{C/T}{\\sinh(C/T)}\\right)^2 + D\\left(\\frac{E/T}{\\cosh(E/T)}\\right)^2$$\n",
"and generalized Plank-Einstein term is given by\n",
"$$\\frac{c_p^0}{R_u} = - \\frac{a_{k} b_{k}^{2} c_{k} d_{k} \\tau^{2} e^{\\frac{b_{k} \\tau}{T_{c}}}}{T_{c}^{2} \\left(c_{k} + d_{k} e^{\\frac{b_{k} \\tau}{T_{c}}}\\right)^{2}}$$\n",
"\n",
"Two important identities:\n",
"$$ \\left(\\frac{1}{\\cosh x}\\right) = \\frac{2e^{-x}}{1+e^{-2x}}$$\n",
"$$ \\left(\\frac{1}{\\sinh x}\\right) = \\frac{2e^{-x}}{1-e^{-2x}}$$\n",
"which yield\n",
"$$ \\left(\\frac{1}{\\cosh x}\\right)^2 = \\frac{4e^{-2x}}{(1+e^{-2x})^2}$$\n",
"$$ \\left(\\frac{1}{\\sinh x}\\right)^2 = \\frac{4e^{-2x}}{(1-e^{-2x})^2}$$"
"Constant is separated out, and handled separately. sinh part can be expanded as\n",
"$$B\\left(\\frac{C/T}{\\sinh(C/T)}\\right)^2 = \\frac{B(-2C/T)^2\\exp(-2C/T)}{(1-\\exp(-2C/T))^2}$$\n",
"where\n",
"$$a_k = B$$\n",
"$$b_k = -2C$$\n",
"$$c_k = 1$$\n",
"$$d_k = -1$$\n",
"\n",
"cosh part can be expanded as\n",
"$$D\\left(\\frac{E/T}{\\cosh(E/T)}\\right)^2 = \\frac{D(-2E/T)^2\\exp(-2E/T)}{(1+\\exp(-2E/T))^2}$$\n",
"where\n",
"$$a_k = -D$$\n",
"$$b_k = -2E$$\n",
"$$c_k = 1$$\n",
"$$d_k = 1$$"
]
},
{
@@ -102,7 +220,7 @@
"display(diff1.rewrite(exp).subs(B, 2*D).simplify())\n",
"term2 = (2*D/T)**2*exp(2*D/T)/(exp(2*D/T)-1)**2\n",
"display(term2)\n",
"print 'which is the same as the Aly-Lee sinh term with b_k = 2*D'\n",
"print 'which is the same as the Aly-Lee sinh term with a_k = C, b_k = 2*D, c_k = -1, d_k = 1'\n",
"\n",
"B, F, T = symbols('B, F, T', real = True, positive = True)\n",
"# From Aly-Lee:\n",

View File

@@ -14,25 +14,24 @@ Residual Helmholtz Energy Terms:
Term | Helmholtz Energy Contribution
---------- | ------------------------------
ResidualHelmholtzPower | \f$ \alpha_r=\left\lbrace\begin{array}{cc}\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} & l_i=0\\ \displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\delta^{l_i}) & l_i\neq 0\end{array}\right.\f$
ResidualHelmholtzExponential | \f$ \alpha_r=\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\gamma_i\delta^{l_i}) \f$
ResidualHelmholtzGaussian | \f$ \alpha_r=\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\eta_i(\delta-\epsilon_i)^2-\beta_i(\tau-\gamma_i)^2)\f$
ResidualHelmholtzGERG2008Gaussian | \f$ \alpha_r=\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\eta_i(\delta-\epsilon_i)^2-\beta_i(\delta-\gamma_i))\f$
ResidualHelmholtzLemmon2005 | \f$ \alpha_r=\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\delta^{l_i}) \exp(-\tau^{m_i})\f$
ResidualHelmholtzNonAnalytic | \f$ \begin{array}{c}\alpha_r&=&\displaystyle\sum_i n_i \Delta^{b_i}\delta\psi \\ \Delta & = & \theta^2+B_i[(\delta-1)^2]^{a_i}\\ \theta & = & (1-\tau)+A_i[(\delta-1)^2]^{1/(2\beta_i)}\\ \psi & = & \exp(-C_i(\delta-1)^2-D_i(\tau-1)^2) \end{array}\f$
ResidualHelmholtzSAFTAssociating | \f$ \alpha_r = am\left(\ln X-\frac{X}{2}+\frac{1}{2}\right); \f$
ResidualHelmholtzPower | \f$ \alpha^r=\left\lbrace\begin{array}{cc}\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} & l_i=0\\ \displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\delta^{l_i}) & l_i\neq 0\end{array}\right.\f$
ResidualHelmholtzExponential | \f$ \alpha^r=\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\gamma_i\delta^{l_i}) \f$
ResidualHelmholtzGaussian | \f$ \alpha^r=\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\eta_i(\delta-\epsilon_i)^2-\beta_i(\tau-\gamma_i)^2)\f$
ResidualHelmholtzGERG2008Gaussian | \f$ \alpha^r=\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\eta_i(\delta-\epsilon_i)^2-\beta_i(\delta-\gamma_i))\f$
ResidualHelmholtzLemmon2005 | \f$ \alpha^r=\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\delta^{l_i}) \exp(-\tau^{m_i})\f$
ResidualHelmholtzNonAnalytic | \f$ \begin{array}{c}\alpha^r&=&\displaystyle\sum_i n_i \Delta^{b_i}\delta\psi \\ \Delta & = & \theta^2+B_i[(\delta-1)^2]^{a_i}\\ \theta & = & (1-\tau)+A_i[(\delta-1)^2]^{1/(2\beta_i)}\\ \psi & = & \exp(-C_i(\delta-1)^2-D_i(\tau-1)^2) \end{array}\f$
ResidualHelmholtzSAFTAssociating | \f$ \alpha^r = am\left(\ln X-\frac{X}{2}+\frac{1}{2}\right); \f$
Ideal-Gas Helmholtz Energy Terms:
Term | Helmholtz Energy Contribution
---------- | ------------------------------
IdealHelmholtzLead | \f$ \alpha_0 = n_1 + n_2\tau + \ln\delta \f$
IdealHelmholtzEnthalpyEntropyOffset | \f$ \alpha_0 = \displaystyle\frac{\Delta s}{R_u/M}+\frac{\Delta h}{(R_u/M)T}\tau \f$
IdealHelmholtzLogTau | \f$ \alpha_0 = n_1\log\tau \f$
IdealHelmholtzPower | \f$ \alpha_0 = \displaystyle\sum_i n_i\tau^{t_i} \f$
IdealHelmholtzPlanckEinstein | \f$ \alpha_0 = \displaystyle\sum_i n_i\log[1-\exp(-\theta_i\tau)] \f$
IdealHelmholtzPlanckEinstein2 | \f$ \alpha_0 = \displaystyle\sum_i n_i\log[c_i+\exp(\theta_i\tau)] \f$
IdealHelmholtzLead | \f$ \alpha^0 = n_1 + n_2\tau + \ln\delta \f$
IdealHelmholtzEnthalpyEntropyOffset | \f$ \alpha^0 = \displaystyle\frac{\Delta s}{R_u/M}+\frac{\Delta h}{(R_u/M)T}\tau \f$
IdealHelmholtzLogTau | \f$ \alpha^0 = n_1\log\tau \f$
IdealHelmholtzPower | \f$ \alpha^0 = \displaystyle\sum_i n_i\tau^{t_i} \f$
IdealHelmholtzPlanckEinsteinGeneralized | \f$ \alpha^0 = \displaystyle\sum_i n_i\log[c_i+d_i\exp(\theta_i\tau)] \f$
*/
class BaseHelmholtzTerm{
public:
@@ -107,7 +106,7 @@ struct ResidualHelmholtzPowerElement
/*!
Term of the form
\f[
\alpha_r=\left\lbrace\begin{array}{cc}\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} & l_i=0\\ \displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\delta^{l_i}) & l_i\neq 0\end{array}\right.
\alpha^r=\left\lbrace\begin{array}{cc}\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} & l_i=0\\ \displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\delta^{l_i}) & l_i\neq 0\end{array}\right.
\f]
*/
class ResidualHelmholtzPower : public BaseHelmholtzTerm{
@@ -162,7 +161,7 @@ struct ResidualHelmholtzExponentialElement
};
/**
Term of the form
\f[ \alpha_r=\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\gamma_i\delta^{l_i}) \f]
\f[ \alpha^r=\displaystyle\sum_i n_i \delta^{d_i} \tau^{t_i} \exp(-\gamma_i\delta^{l_i}) \f]
*/
class ResidualHelmholtzExponential : public BaseHelmholtzTerm{
@@ -573,7 +572,7 @@ public:
/// The leading term in the EOS used to set the desired reference state
/**
\f[
\alpha_0 = \log(\delta)+a_1+a_2\tau
\alpha^0 = \log(\delta)+a_1+a_2\tau
\f]
*/
class IdealHelmholtzLead : public BaseHelmholtzTerm{
@@ -632,7 +631,7 @@ public:
/// The term in the EOS used to shift the reference state of the fluid
/**
\f[
\alpha_0 = a_1+a_2\tau
\alpha^0 = a_1+a_2\tau
\f]
*/
class IdealHelmholtzEnthalpyEntropyOffset : public BaseHelmholtzTerm{
@@ -682,7 +681,7 @@ public:
/**
\f[
\alpha_0 = a_1\ln\tau
\alpha^0 = a_1\ln\tau
\f]
*/
class IdealHelmholtzLogTau : public BaseHelmholtzTerm
@@ -735,7 +734,7 @@ public:
/**
\f[
\alpha_0 = \displaystyle\sum_i n_i\tau^{t_i}
\alpha^0 = \displaystyle\sum_i n_i\tau^{t_i}
\f]
*/
class IdealHelmholtzPower : public BaseHelmholtzTerm{
@@ -793,33 +792,42 @@ public:
/**
\f[
\alpha_0 = \displaystyle\sum_i n_i\log[1-\exp(-\theta_i\tau)]
\alpha^0 = \displaystyle\sum_i n_i\log[c_i+d_i\exp(\theta_i\tau)]
\f]
*/
class IdealHelmholtzPlanckEinstein : public BaseHelmholtzTerm{
class IdealHelmholtzPlanckEinsteinGeneralized : public BaseHelmholtzTerm{
private:
std::vector<long double> n,theta; // Use these variables internally
std::vector<long double> n,theta,c,d; // Use these variables internally
std::size_t N;
bool enabled;
public:
IdealHelmholtzPlanckEinstein(){N = 0; enabled = false;}
IdealHelmholtzPlanckEinsteinGeneralized(){N = 0; enabled = false;}
// Constructor with std::vector instances
IdealHelmholtzPlanckEinstein(std::vector<long double> n, std::vector<long double> theta)
:n(n), theta(theta)
IdealHelmholtzPlanckEinsteinGeneralized(std::vector<long double> n, std::vector<long double> theta, std::vector<long double> c, std::vector<long double> d)
:n(n), theta(theta), c(c), d(d)
{
N = n.size();
enabled = true;
};
//Destructor
~IdealHelmholtzPlanckEinstein(){};
// Destructor
~IdealHelmholtzPlanckEinsteinGeneralized(){};
// Extend the vectors to allow for multiple instances feeding values to this function
void extend(std::vector<long double> n, std::vector<long double> theta, std::vector<long double> c, std::vector<long double> d)
{
this->n.insert(this->n.end(), n.begin(), n.end());
this->theta.insert(this->theta.end(), theta.begin(), theta.end());
this->c.insert(this->c.end(), c.begin(), c.end());
this->d.insert(this->d.end(), d.begin(), d.end());
}
bool is_enabled(){return enabled;};
void to_json(rapidjson::Value &el, rapidjson::Document &doc)
{
el.AddMember("type","IdealHelmholtzPlanckEinstein",doc.GetAllocator());
el.AddMember("type","IdealHelmholtzPlanckEinsteinGeneralized",doc.GetAllocator());
cpjson::set_long_double_array("n",n,el,doc);
cpjson::set_long_double_array("theta",theta,el,doc);
};
@@ -828,81 +836,22 @@ public:
long double base(const long double &tau, const long double &delta) throw(){
if (!enabled){return 0.0;}
long double s=0; for (std::size_t i=0; i < N; ++i){
s += n[i]*log(1.0-exp(-theta[i]*tau));
} return s;
s += n[i]*log(c[i]+d[i]*exp(theta[i]*tau));
}
return s;
};
long double dTau(const long double &tau, const long double &delta) throw(){
if (!enabled){return 0.0;}
long double s=0; for (std::size_t i=0; i < N; ++i){s += n[i]*theta[i]*(1.0/(1.0-exp(-theta[i]*tau))-1.0);} return s;
long double s=0; for (std::size_t i=0; i < N; ++i){s += n[i]*theta[i]*d[i]*exp(theta[i]*tau)/(c[i]+d[i]*exp(theta[i]*tau));}
return s;
};
long double dTau2(const long double &tau, const long double &delta) throw(){
if (!enabled){return 0.0;}
long double s=0; for (std::size_t i=0; i < N; ++i){s -= n[i]*pow(theta[i],2)*exp(theta[i]*tau)/pow(1.0-exp(theta[i]*tau),2);} return s;
long double s=0; for (std::size_t i=0; i < N; ++i){s += n[i]*pow(theta[i],2)*c[i]*d[i]*exp(theta[i]*tau)/pow(c[i]+d[i]*exp(theta[i]*tau),2);} return s;
};
long double dTau3(const long double &tau, const long double &delta) throw(){
if (!enabled){return 0.0;}
long double s=0; for (std::size_t i=0; i < N; ++i){s += n[i]*pow(theta[i],2)*theta[i]*exp(theta[i]*tau)*(exp(theta[i]*tau)+1)/pow(exp(theta[i]*tau)-1,3);} return s;
};
long double dDelta(const long double &tau, const long double &delta) throw(){return 0.0;};
long double dDelta2(const long double &tau, const long double &delta) throw(){return 0.0;};
long double dDelta2_dTau(const long double &tau, const long double &delta) throw(){return 0.0;};
long double dDelta_dTau(const long double &tau, const long double &delta) throw(){return 0.0;};
long double dDelta_dTau2(const long double &tau, const long double &delta) throw(){return 0.0;};
long double dDelta3(const long double &tau, const long double &delta) throw(){return 0;};
};
/**
\f[
\alpha_0 = \displaystyle\sum_i n_i\log[c_i+\exp(\theta_i\tau)]
\f]
*/
class IdealHelmholtzPlanckEinstein2 : public BaseHelmholtzTerm{
private:
std::vector<long double> n,theta,c; // Use these variables internally
std::size_t N;
bool enabled;
public:
IdealHelmholtzPlanckEinstein2(){N = 0; enabled = false;}
// Constructor with std::vector instances
IdealHelmholtzPlanckEinstein2(const std::vector<long double> &n,
const std::vector<long double> &theta,
const std::vector<long double> &c)
:n(n), theta(theta), c(c)
{
N = n.size();
enabled = true;
};
//Destructor
~IdealHelmholtzPlanckEinstein2(){};
bool is_enabled(){return enabled;};
void to_json(rapidjson::Value &el, rapidjson::Document &doc)
{
el.AddMember("type","IdealHelmholtzPlanckEinstein2",doc.GetAllocator());
cpjson::set_long_double_array("n",n,el,doc);
cpjson::set_long_double_array("theta",theta,el,doc);
cpjson::set_long_double_array("c",c,el,doc);
};
// Term and its derivatives
long double base(const long double &tau, const long double &delta) throw(){
if (!enabled){return 0.0;}
long double s=0; for (std::size_t i=0; i < N; ++i){s += n[i]*log(c[i]+exp(theta[i]*tau));} return s;
};
long double dTau(const long double &tau, const long double &delta) throw(){
if (!enabled){return 0.0;}
long double s=0; for (std::size_t i=0; i < N; ++i){s += n[i]*theta[i]*exp(tau*theta[i])/(c[i]+exp(theta[i]*tau));} return s;
};
long double dTau2(const long double &tau, const long double &delta) throw(){
if (!enabled){return 0.0;}
long double s=0; for (std::size_t i=0; i < N; ++i){s += n[i]*pow(theta[i],2)*c[i]*exp(tau*theta[i])/pow(c[i]+exp(tau*theta[i]),2);} return s;
};
long double dTau3(const long double &tau, const long double &delta) throw(){
if (!enabled){return 0.0;}
long double s=0; for (std::size_t i=0; i < N; ++i){s += n[i]*pow(theta[i],2)*c[i]*(-theta[i])*exp(theta[i]*tau)*(exp(theta[i]*tau)-c[i])/pow(exp(theta[i]*tau)+c[i],3);} return s;
long double s=0; for (std::size_t i=0; i < N; ++i){s += n[i]*pow(theta[i],3)*c[i]*d[i]*(c[i]-d[i]*exp(theta[i]*tau))*exp(theta[i]*tau)/pow(c[i]+d[i]*exp(theta[i]*tau),3);} return s;
};
long double dDelta(const long double &tau, const long double &delta) throw(){return 0.0;};
long double dDelta2(const long double &tau, const long double &delta) throw(){return 0.0;};
@@ -986,6 +935,13 @@ public:
N = c.size();
};
void extend(const std::vector<long double> &c, const std::vector<long double> &t)
{
this->c.insert(this->c.end(), c.begin(), c.end());
this->t.insert(this->t.end(), t.begin(), t.end());
N += c.size();
}
/// Destructor
~IdealHelmholtzCP0PolyT(){};
@@ -1147,100 +1103,89 @@ public:
IdealHelmholtzEnthalpyEntropyOffset EnthalpyEntropyOffset;
IdealHelmholtzLogTau LogTau;
IdealHelmholtzPower Power;
IdealHelmholtzPlanckEinstein PlanckEinstein;
IdealHelmholtzPlanckEinstein2 PlanckEinstein2;
IdealHelmholtzPlanckEinsteinGeneralized PlanckEinstein;
IdealHelmholtzCP0Constant CP0Constant;
IdealHelmholtzCP0PolyT CP0PolyT;
IdealHelmholtzCP0AlyLee CP0AlyLee;
long double base(const long double &tau, const long double &delta)
{
return (Lead.base(tau, delta) + EnthalpyEntropyOffset.base(tau, delta)
+ LogTau.base(tau, delta) + Power.base(tau, delta)
+ PlanckEinstein.base(tau, delta) + PlanckEinstein2.base(tau, delta)
+ PlanckEinstein.base(tau, delta)
+ CP0Constant.base(tau, delta) + CP0PolyT.base(tau, delta)
+ CP0AlyLee.base(tau, delta)
);
};
long double dDelta(const long double &tau, const long double &delta)
{
return (Lead.dDelta(tau, delta) + EnthalpyEntropyOffset.dDelta(tau, delta)
+ LogTau.dDelta(tau, delta) + Power.dDelta(tau, delta)
+ PlanckEinstein.dDelta(tau, delta) + PlanckEinstein2.dDelta(tau, delta)
+ PlanckEinstein.dDelta(tau, delta)
+ CP0Constant.dDelta(tau, delta) + CP0PolyT.dDelta(tau, delta)
+ CP0AlyLee.dDelta(tau, delta)
);
};
long double dTau(const long double &tau, const long double &delta)
{
return (Lead.dTau(tau, delta) + EnthalpyEntropyOffset.dTau(tau, delta)
+ LogTau.dTau(tau, delta) + Power.dTau(tau, delta)
+ PlanckEinstein.dTau(tau, delta) + PlanckEinstein2.dTau(tau, delta)
+ PlanckEinstein.dTau(tau, delta)
+ CP0Constant.dTau(tau, delta) + CP0PolyT.dTau(tau, delta)
+ CP0AlyLee.dTau(tau, delta)
);
};
long double dDelta2(const long double &tau, const long double &delta)
{
return (Lead.dDelta2(tau, delta) + EnthalpyEntropyOffset.dDelta2(tau, delta)
+ LogTau.dDelta2(tau, delta) + Power.dDelta2(tau, delta)
+ PlanckEinstein.dDelta2(tau, delta) + PlanckEinstein2.dDelta2(tau, delta)
+ PlanckEinstein.dDelta2(tau, delta)
+ CP0Constant.dDelta2(tau, delta) + CP0PolyT.dDelta2(tau, delta)
+ CP0AlyLee.dDelta2(tau, delta)
);
};
long double dDelta_dTau(const long double &tau, const long double &delta)
{
return (Lead.dDelta_dTau(tau, delta) + EnthalpyEntropyOffset.dDelta_dTau(tau, delta)
+ LogTau.dDelta_dTau(tau, delta) + Power.dDelta_dTau(tau, delta)
+ PlanckEinstein.dDelta_dTau(tau, delta) + PlanckEinstein2.dDelta_dTau(tau, delta)
+ PlanckEinstein.dDelta_dTau(tau, delta)
+ CP0Constant.dDelta_dTau(tau, delta) + CP0PolyT.dDelta_dTau(tau, delta)
+ CP0AlyLee.dDelta_dTau(tau, delta)
);
};
long double dTau2(const long double &tau, const long double &delta)
{
return (Lead.dTau2(tau, delta) + EnthalpyEntropyOffset.dTau2(tau, delta)
+ LogTau.dTau2(tau, delta) + Power.dTau2(tau, delta)
+ PlanckEinstein.dTau2(tau, delta) + PlanckEinstein2.dTau2(tau, delta)
+ PlanckEinstein.dTau2(tau, delta)
+ CP0Constant.dTau2(tau, delta) + CP0PolyT.dTau2(tau, delta)
+ CP0AlyLee.dTau2(tau, delta));
);
};
long double dDelta3(const long double &tau, const long double &delta)
{
return (Lead.dDelta3(tau, delta) + EnthalpyEntropyOffset.dDelta3(tau, delta)
+ LogTau.dDelta3(tau, delta) + Power.dDelta3(tau, delta)
+ PlanckEinstein.dDelta3(tau, delta) + PlanckEinstein2.dDelta3(tau, delta)
+ PlanckEinstein.dDelta3(tau, delta)
+ CP0Constant.dDelta3(tau, delta) + CP0PolyT.dDelta3(tau, delta)
+ CP0AlyLee.dDelta3(tau, delta)
);
};
long double dDelta2_dTau(const long double &tau, const long double &delta)
{
return (Lead.dDelta2_dTau(tau, delta) + EnthalpyEntropyOffset.dDelta2_dTau(tau, delta)
+ LogTau.dDelta2_dTau(tau, delta) + Power.dDelta2_dTau(tau, delta)
+ PlanckEinstein.dDelta2_dTau(tau, delta) + PlanckEinstein2.dDelta2_dTau(tau, delta)
+ PlanckEinstein.dDelta2_dTau(tau, delta)
+ CP0Constant.dDelta2_dTau(tau, delta) + CP0PolyT.dDelta2_dTau(tau, delta)
+ CP0AlyLee.dDelta2_dTau(tau, delta)
);
};
long double dDelta_dTau2(const long double &tau, const long double &delta)
{
return (Lead.dDelta_dTau2(tau, delta) + EnthalpyEntropyOffset.dDelta_dTau2(tau, delta)
+ LogTau.dDelta_dTau2(tau, delta) + Power.dDelta_dTau2(tau, delta)
+ PlanckEinstein.dDelta_dTau2(tau, delta) + PlanckEinstein2.dDelta_dTau2(tau, delta)
+ PlanckEinstein.dDelta_dTau2(tau, delta)
+ CP0Constant.dDelta_dTau2(tau, delta) + CP0PolyT.dDelta_dTau2(tau, delta)
+ CP0AlyLee.dDelta_dTau2(tau, delta)
);
};
long double dTau3(const long double &tau, const long double &delta)
{
return (Lead.dTau3(tau, delta) + EnthalpyEntropyOffset.dTau3(tau, delta)
+ LogTau.dTau3(tau, delta) + Power.dTau3(tau, delta)
+ PlanckEinstein.dTau3(tau, delta) + PlanckEinstein2.dTau3(tau, delta)
+ PlanckEinstein.dTau3(tau, delta)
+ CP0Constant.dTau3(tau, delta) + CP0PolyT.dTau3(tau, delta)
+ CP0AlyLee.dTau3(tau, delta)
);
};
};

View File

@@ -143,20 +143,36 @@ protected:
long double a = cpjson::get_double(contribution,"a");
EOS.alpha0.LogTau = IdealHelmholtzLogTau(a);
}
else if (!type.compare("IdealGasHelmholtzPlanckEinsteinGeneralized"))
{
// Retrieve the values
std::vector<long double> n = cpjson::get_long_double_array(contribution["n"]);
std::vector<long double> t = cpjson::get_long_double_array(contribution["t"]);
std::vector<long double> c = cpjson::get_long_double_array(contribution["c"]);
std::vector<long double> d = cpjson::get_long_double_array(contribution["d"]);
if (EOS.alpha0.PlanckEinstein.is_enabled() == true){
EOS.alpha0.PlanckEinstein.extend(n, t, c, d);
}
else{
EOS.alpha0.PlanckEinstein = IdealHelmholtzPlanckEinsteinGeneralized(n, t, c, d);
}
}
else if (!type.compare("IdealGasHelmholtzPlanckEinstein"))
{
if (EOS.alpha0.PlanckEinstein.is_enabled() == true){throw ValueError("Cannot add ");}
// Retrieve the values
std::vector<long double> n = cpjson::get_long_double_array(contribution["n"]);
std::vector<long double> t = cpjson::get_long_double_array(contribution["t"]);
EOS.alpha0.PlanckEinstein = IdealHelmholtzPlanckEinstein(n, t);
}
else if (!type.compare("IdealGasHelmholtzPlanckEinstein2"))
{
if (EOS.alpha0.PlanckEinstein2.is_enabled() == true){throw ValueError("Cannot add");}
std::vector<long double> n = cpjson::get_long_double_array(contribution["n"]);
std::vector<long double> t = cpjson::get_long_double_array(contribution["t"]);
std::vector<long double> c = cpjson::get_long_double_array(contribution["c"]);
EOS.alpha0.PlanckEinstein2 = IdealHelmholtzPlanckEinstein2(n, t, c);
// Flip the sign of theta
for (std::size_t i = 0; i < t.size(); ++i){ t[i] *= -1;}
std::vector<long double> c(n.size(), 1);
std::vector<long double> d(c.size(), -1);
if (EOS.alpha0.PlanckEinstein.is_enabled() == true){
EOS.alpha0.PlanckEinstein.extend(n, t, c, d);
}
else{
EOS.alpha0.PlanckEinstein = IdealHelmholtzPlanckEinsteinGeneralized(n, t, c, d);
}
}
else if (!type.compare("IdealGasHelmholtzCP0Constant"))
{
@@ -177,11 +193,45 @@ protected:
}
else if (!type.compare("IdealGasHelmholtzCP0AlyLee"))
{
if (EOS.alpha0.CP0AlyLee.is_enabled() == true){std::cout << "Cannot add IdealGasHelmholtzCP0AlyLee\n";}
std::vector<long double> c = cpjson::get_long_double_array(contribution["c"]);
std::vector<long double> constants = cpjson::get_long_double_array(contribution["c"]);
long double Tc = cpjson::get_double(contribution, "Tc");
long double T0 = cpjson::get_double(contribution, "T0");
EOS.alpha0.CP0AlyLee = IdealHelmholtzCP0AlyLee(c, Tc, T0);
// Take the constant term if nonzero and set it as a polyT term
if (fabs(constants[0]) > 1e-14){
std::vector<long double> c(1,constants[0]), t(1,0);
if (EOS.alpha0.CP0PolyT.is_enabled() == true){
EOS.alpha0.CP0PolyT.extend(c,t);
}
else{
EOS.alpha0.CP0PolyT = IdealHelmholtzCP0PolyT(c, t, Tc, T0);
}
}
std::vector<long double> n, c, d, t;
if (fabs(constants[1]) > 1e-14){
// sinh term can be converted by setting a_k = C, b_k = 2*D, c_k = -1, d_k = 1
n.push_back(constants[1]);
t.push_back(-2*constants[2]/Tc);
c.push_back(1);
d.push_back(-1);
}
if (fabs(constants[3]) > 1e-14){
// cosh term can be converted by setting a_k = C, b_k = 2*D, c_k = 1, d_k = 1
n.push_back(-constants[3]);
t.push_back(-2*constants[4]/Tc);
c.push_back(1);
d.push_back(1);
}
if (EOS.alpha0.PlanckEinstein.is_enabled() == true){
EOS.alpha0.PlanckEinstein.extend(n, t, c, d);
}
else{
EOS.alpha0.PlanckEinstein = IdealHelmholtzPlanckEinsteinGeneralized(n, t, c, d);
}
}
else if (!type.compare("IdealGasHelmholtzEnthalpyEntropyOffset"))
{

View File

@@ -1715,42 +1715,50 @@ long double HelmholtzEOSMixtureBackend::calc_alphar_deriv_nocache(const int nTau
long double HelmholtzEOSMixtureBackend::calc_alpha0_deriv_nocache(const int nTau, const int nDelta, const std::vector<long double> &mole_fractions,
const long double &tau, const long double &delta, const long double &Tr, const long double &rhor)
{
long double val;
if (is_pure_or_pseudopure)
{
if (nTau == 0 && nDelta == 0){
return components[0]->pEOS->base0(tau, delta);
val = components[0]->pEOS->base0(tau, delta);
}
else if (nTau == 0 && nDelta == 1){
return components[0]->pEOS->dalpha0_dDelta(tau, delta);
val = components[0]->pEOS->dalpha0_dDelta(tau, delta);
}
else if (nTau == 1 && nDelta == 0){
return components[0]->pEOS->dalpha0_dTau(tau, delta);
val = components[0]->pEOS->dalpha0_dTau(tau, delta);
}
else if (nTau == 0 && nDelta == 2){
return components[0]->pEOS->d2alpha0_dDelta2(tau, delta);
val = components[0]->pEOS->d2alpha0_dDelta2(tau, delta);
}
else if (nTau == 1 && nDelta == 1){
return components[0]->pEOS->d2alpha0_dDelta_dTau(tau, delta);
val = components[0]->pEOS->d2alpha0_dDelta_dTau(tau, delta);
}
else if (nTau == 2 && nDelta == 0){
return components[0]->pEOS->d2alpha0_dTau2(tau, delta);
val = components[0]->pEOS->d2alpha0_dTau2(tau, delta);
}
else if (nTau == 0 && nDelta == 3){
return components[0]->pEOS->d3alpha0_dDelta3(tau, delta);
val = components[0]->pEOS->d3alpha0_dDelta3(tau, delta);
}
else if (nTau == 1 && nDelta == 2){
return components[0]->pEOS->d3alpha0_dDelta2_dTau(tau, delta);
val = components[0]->pEOS->d3alpha0_dDelta2_dTau(tau, delta);
}
else if (nTau == 2 && nDelta == 1){
return components[0]->pEOS->d3alpha0_dDelta_dTau2(tau, delta);
val = components[0]->pEOS->d3alpha0_dDelta_dTau2(tau, delta);
}
else if (nTau == 3 && nDelta == 0){
return components[0]->pEOS->d3alpha0_dTau3(tau, delta);
val = components[0]->pEOS->d3alpha0_dTau3(tau, delta);
}
else
{
throw ValueError();
}
if (!ValidNumber(val)){
calc_alpha0_deriv_nocache(nTau,nDelta,mole_fractions,tau,delta,Tr,rhor);
throw ValueError(format("calc_alpha0_deriv_nocache returned invalid number with inputs nTau: %d, nDelta: %d", nTau, nDelta));
}
else{
return val;
}
}
else{
// See Table B5, GERG 2008 from Kunz Wagner, JCED, 2012

View File

@@ -1691,10 +1691,11 @@ long double IdealHelmholtzCP0AlyLee::dTau(const long double &tau, const long dou
}
long double IdealHelmholtzCP0AlyLee::anti_deriv_cp0_tau2(const long double &tau)
{
return -c[0]/tau + 2*c[1]*c[2]/(Tc*(exp(-(2*c[2]*tau)/Tc)-1)) + 2*c[3]*c[4]/(Tc*(exp(-(2*c[4]*tau)/Tc)+1));
return -c[0]/tau + 2*c[1]*c[2]/Tc/(exp(-2*c[2]*tau/Tc)-1) - 2*c[3]*c[4]/Tc*(exp(2*c[4]*tau/Tc)+1);
}
long double IdealHelmholtzCP0AlyLee::anti_deriv_cp0_tau(const long double &tau)
{
long double lnarg =-1 + exp(-2*c[2]*tau/Tc);
long double term1 = c[0]*log(tau);
long double term2 = 2*c[1]*c[2]*tau/(-Tc + Tc*exp(-2*c[2]*tau/Tc)) + c[1]*log(-1 + exp(-2*c[2]*tau/Tc)) + 2*c[1]*c[2]*tau/Tc;
long double term3 = -c[3]*(Tc*exp(2*c[4]*tau/Tc)*log(exp(2*c[4]*tau/Tc) + 1) + Tc*log(exp(2*c[4]*tau/Tc) + 1) - 2*c[4]*tau*exp(2*c[4]*tau/Tc))/(Tc*(exp(2*c[4]*tau/Tc) + 1));
@@ -1716,7 +1717,6 @@ long double IdealHelmholtzCP0AlyLee::dTau3(const long double &tau, const long do
/*
IdealHelmholtzEnthalpyEntropyOffset EnthalpyEntropyOffset;
IdealHelmholtzCP0AlyLee CP0AlyLee;
*/
#ifdef ENABLE_CATCH
@@ -1728,8 +1728,8 @@ class HelmholtzConsistencyFixture
public:
long double numerical, analytic;
std::tr1::shared_ptr<CoolProp::BaseHelmholtzTerm> Lead, LogTau, IGPower, PlanckEinstein, PlanckEinstein2,
CP0Constant, CP0PolyT, CP0AlyLee, Gaussian, Lemmon2005, Power, SAFT, NonAnalytic, Exponential, GERG2008;
std::tr1::shared_ptr<CoolProp::BaseHelmholtzTerm> Lead, LogTau, IGPower, PlanckEinstein,
CP0Constant, CP0PolyT, Gaussian, Lemmon2005, Power, SAFT, NonAnalytic, Exponential, GERG2008;
HelmholtzConsistencyFixture(){
Lead.reset(new CoolProp::IdealHelmholtzLead(1,3));
@@ -1739,12 +1739,8 @@ public:
IGPower.reset(new CoolProp::IdealHelmholtzPower(n,t));
}
{
std::vector<long double> n(4,0), t(4,1); n[0] = 0.1; n[2] = 0.5; t[1] = 1; t[2] = 2; t[3] = 2;
PlanckEinstein.reset(new CoolProp::IdealHelmholtzPlanckEinstein(n, t));
}
{
std::vector<long double> n(4,0), t(4,1), c(4,1); n[0] = -0.1; n[2] = 0.1; t[1] = -1; t[2] = -2; t[3] = 2;
PlanckEinstein2.reset(new CoolProp::IdealHelmholtzPlanckEinstein2(n, t, c));
std::vector<long double> n(4,0), t(4,1), c(4,1), d(4,1); n[0] = 0.1; n[2] = 0.5; t[1] = 1; t[2] = 2; t[3] = 2;
PlanckEinstein.reset(new CoolProp::IdealHelmholtzPlanckEinsteinGeneralized(n, t, c, d));
}
{
long double T0 = 273.15, Tc = 345.857, c = 1.0578, t = 0.33;
@@ -1752,11 +1748,6 @@ public:
std::vector<long double>(1,t),
Tc, T0));
}
{
long double T0 = 518.109977174843, Tc = 645.78, c[] = {56.37158920013201, 118.0111016069331, 1792.1, 82.5909330141469, 786.8};
CP0AlyLee.reset(new CoolProp::IdealHelmholtzCP0AlyLee(std::vector<long double>(c, c+sizeof(c)/sizeof(c[0])),
Tc, T0));
}
CP0Constant.reset(new CoolProp::IdealHelmholtzCP0Constant(4/8.314472,300,250));
{
long double beta[] = {1.24, 0.821, 15.45, 2.21, 437, 0.743},
@@ -1871,10 +1862,8 @@ public:
else if (!t.compare("LogTau")){return LogTau;}
else if (!t.compare("IGPower")){return IGPower;}
else if (!t.compare("PlanckEinstein")){return PlanckEinstein;}
else if (!t.compare("PlanckEinstein2")){return PlanckEinstein2;}
else if (!t.compare("CP0Constant")){return CP0Constant;}
else if (!t.compare("CP0PolyT")){return CP0PolyT;}
else if (!t.compare("CP0AlyLee")){return CP0AlyLee;}
else if (!t.compare("Gaussian")){return Gaussian;}
else if (!t.compare("Lemmon2005")){return Lemmon2005;}

View File

@@ -154,8 +154,8 @@ int main()
}
if (1)
{
double rrr0 = PropsSI("O","T",350,"Q",0,"REFPROP::Water");
double rrr2 = PropsSI("O","T",350,"Q",0,"Water");
double rrr0 = PropsSI("C","T",350,"D",1e-13,"REFPROP::MDM");
double rrr2 = PropsSI("C","T",350,"D",1e-13,"MDM");
double rrr =0 ;
}
if (1)
@@ -171,7 +171,8 @@ int main()
}
if (1)
{
double h1 = PropsSI("H","T",273.15,"Q",0,"CO2");
double h1 = PropsSI("S","P",101325,"Q",0,"n-Pentane");
std::string er = get_global_param_string("errstring");
set_reference_stateS("n-Propane","NBP");
double h2 = PropsSI("H","P",101325,"Q",0,"n-Propane");

View File

@@ -72,16 +72,6 @@ cpdef ndarray_or_iterable(object input):
# if retval < 0:
# raise ValueError('Unable to set reference state')
#
# cpdef add_REFPROP_fluid(string_like FluidName):
# """
# Add a REFPROP fluid to CoolProp internal structure
#
# example::
#
# add_REFPROP_fluid("REFPROP-PROPANE")
#
# """
# _add_REFPROP_fluid(FluidName)
#
# cpdef long get_Fluid_index(string_like Fluid):
# """
@@ -116,9 +106,9 @@ cpdef ndarray_or_iterable(object input):
cpdef get_global_param_string(string param):
return _get_global_param_string(param)
#
# cpdef get_fluid_param_string(string_like fluid, string_like param):
# return _get_fluid_param_string(fluid, param)
cpdef get_fluid_param_string(string_like fluid, string_like param):
return _get_fluid_param_string(fluid, param)
#
# cpdef __Props_err1(in1,in2,errstr):
# if not len(errstr) == 0:
@@ -155,53 +145,53 @@ cpdef PropsSI(in1, in2, in3, in4, in5, in6, in7 = None):
else:
return _PropsSI(in1, in2, in3, in4, in5, in6, in7)
# cpdef list FluidsList():
# """
# Return a list of strings of all fluid names
#
# Returns
# -------
# FluidsList : list of strings of fluid names
# All the fluids that are included in CoolProp
#
# Notes
# -----
#
# Here is an example::
#
# In [0]: from CoolProp.CoolProp import FluidsList
#
# In [1]: FluidsList()
#
# """
# return _get_global_param_string("FluidsList").split(',')
cpdef list FluidsList():
"""
Return a list of strings of all fluid names
Returns
-------
FluidsList : list of strings of fluid names
All the fluids that are included in CoolProp
Notes
-----
Here is an example::
In [0]: from CoolProp.CoolProp import FluidsList
In [1]: FluidsList()
"""
return _get_global_param_string("FluidsList").split(',')
# cpdef get_aliases(str Fluid):
# """
# Return a comma separated string of aliases for the given fluid
# """
# cdef bytes _Fluid = Fluid.encode('ascii')
# return [F.encode('ascii') for F in (_get_fluid_param_string(_Fluid,'aliases').encode('ascii')).decode('ascii').split(',')]
cpdef get_aliases(str Fluid):
"""
Return a comma separated string of aliases for the given fluid
"""
cdef bytes _Fluid = Fluid.encode('ascii')
return [F.encode('ascii') for F in (_get_fluid_param_string(_Fluid,'aliases').encode('ascii')).decode('ascii').split(',')]
#
# cpdef string get_REFPROPname(str Fluid):
# """
# Return the REFPROP compatible name for the fluid (only useful on windows)
#
# Some fluids do not use the REFPROP name. For instance,
# ammonia is R717, and propane is R290. You can still can still call CoolProp
# using the name ammonia or R717, but REFPROP requires that you use a limited
# subset of names. Therefore, this function that returns the REFPROP compatible
# name. To then use this to call REFPROP, you would do something like::
#
# In [0]: from CoolProp.CoolProp import get_REFPROPname, Props
#
# In [1]: Fluid = 'REFPROP-' + get_REFPROPname('R290')
#
# In [2]: Props('D', 'T', 300, 'P', 300, Fluid)
# """
# cdef bytes _Fluid = Fluid.encode('ascii')
# return _get_fluid_param_string(_Fluid,'REFPROP_name')
cpdef string get_REFPROPname(str Fluid):
"""
Return the REFPROP compatible name for the fluid
Some fluids do not use the REFPROP name. For instance,
ammonia is R717, and propane is R290. You can still can still call CoolProp
using the name ammonia or R717, but REFPROP requires that you use a limited
subset of names. Therefore, this function that returns the REFPROP compatible
name. To then use this to call REFPROP, you would do something like::
In [0]: from CoolProp.CoolProp import get_REFPROPname, PropsSI
In [1]: get_REFPROPname('R290')
In [2]: PropsSI('D', 'T', 300, 'P', 300, Fluid)
"""
cdef bytes _Fluid = Fluid.encode('ascii')
return _get_fluid_param_string(_Fluid,'REFPROP_name')
# cpdef string get_BibTeXKey(str Fluid, str key):
# """
@@ -230,31 +220,31 @@ cpdef PropsSI(in1, in2, in3, in4, in5, in6, in7 = None):
# Return the current error string
# """
# return _get_global_param_string("errstring")
# cpdef set_debug_level(int level):
# """
# Set the current debug level as integer in the range [0,10]
#
# Parameters
# ----------
# level : int
# If level is 0, no output will be written to screen, if >0,
# some output will be written to screen. The larger level is,
# the more verbose the output will be
# """
# _set_debug_level(level)
cpdef set_debug_level(int level):
"""
Set the current debug level as integer in the range [0,10]
Parameters
----------
level : int
If level is 0, no output will be written to screen, if >0,
some output will be written to screen. The larger level is,
the more verbose the output will be
"""
_set_debug_level(level)
# cpdef get_debug_level():
# """
# Return the current debug level as integer
#
# Returns
# -------
# level : int
# If level is 0, no output will be written to screen, if >0,
# some output will be written to screen. The larger level is,
# the more verbose the output will be
# """
# return _get_debug_level()
cpdef int get_debug_level():
"""
Return the current debug level as integer
Returns
-------
level : int
If level is 0, no output will be written to screen, if >0,
some output will be written to screen. The larger level is,
the more verbose the output will be
"""
return _get_debug_level()
#
# cpdef bint IsFluidType(string Fluid, string Type):
# """
@@ -286,8 +276,8 @@ cdef dict paras = {iDmass : 'D',
iCpmass : 'C',
# iC0 : 'C0',
iCvmass : 'O',
# iV : 'V',
# iL : 'L',
iV : 'V',
iL : 'L',
iSmass : 'S',
iUmass : 'U',
# iDpdT : 'dpdT'
@@ -307,16 +297,16 @@ cdef class State:
to calculate the enthalpy and pressure. Since the Equations of State are
all explicit in temperature and density, each time you call something like::
h = Props('H','T',T','P',P,Fluid)
s = Props('S','T',T','P',P,Fluid)
h = PropsSI('H','T',T','P',P,Fluid)
s = PropsSI('S','T',T','P',P,Fluid)
the solver is used to carry out the T-P flash calculation. And if you wanted
entropy as well you could either intermediately calculate ``T``, ``rho`` and then use
``T``, ``rho`` in the EOS in a manner like::
rho = Props('D','T',T','P',P,Fluid)
h = Props('H','T',T','D',rho,Fluid)
s = Props('S','T',T','D',rho,Fluid)
rho = PropsSI('D','T',T','P',P,Fluid)
h = PropsSI('H','T',T','D',rho,Fluid)
s = PropsSI('S','T',T','D',rho,Fluid)
Instead in this class all that is handled internally. So the call to update
sets the internal variables in the most computationally efficient way possible
@@ -332,7 +322,7 @@ cdef class State:
phase : string
DEPRECATED : this input is ignored
backend : string
The CoolProp backend that should be used
The CoolProp backend that should be used, one of "HEOS" (default), "REFPROP", "INCOMP", "BRINE", etc.
"""
cdef string _Fluid = Fluid