{
 "metadata": {
  "name": "",
  "signature": "sha256:fa5ab579b67d03b137ac3fa2a3c15807169799511a42efb187eb3ac12b750898"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "heading",
     "level": 1,
     "metadata": {},
     "source": [
      "Forming the Phase Envelope"
     ]
    },
    {
     "cell_type": "heading",
     "level": 2,
     "metadata": {},
     "source": [
      "Overview"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The analysis in this section follows the methodologies proposed in the GERG 2004 monograph published in 2007\n",
      "\n",
      "System of Equations"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Our residual vector $\\mathbf{F}$ is equal to \n",
      "\n",
      "$$F_i = \\ln\\phi(T,p,\\mathbf{y})-\\ln \\phi(T,p,\\mathbf{x})+\\ln K_i=0,  i=1,2,3... N$$\n",
      "\n",
      "$$F_{N+1} = \\sum_{i=1}^{N}(y_i-x_i)=0$$\n",
      "    \n",
      "where\n",
      "\n",
      "$$x_i = \\frac{z_i}{1-\\beta+\\beta K_i}$$\n",
      "    \n",
      "and \n",
      "\n",
      "$$y_i = \\frac{K_iz_i}{1-\\beta+\\beta K_i}$$\n",
      "    \n",
      "DO NOT NORMALIZE $x$ and $y$ !!!!\n",
      "\n",
      "$$F_{N+2} = X_s - S = 0$$\n",
      "    \n",
      "and the system to be solved is equal to\n",
      "\n",
      "$$\\mathbf{J}\\mathbf{\\Delta X}= -\\mathbf{F}$$\n",
      "    "
     ]
    },
    {
     "cell_type": "heading",
     "level": 2,
     "metadata": {},
     "source": [
      "Building the Jacobian matrix"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "This is the trickiest part of this method.  There are a lot of derivatives to implement, and we want to implement all of them analytically.\n",
      "\n",
      "$$\\frac{\\partial F_i}{\\partial \\ln T} = T\\left[ \\left(\\frac{\\partial \\ln \\phi_i}{\\partial T}\\right)''_{p,\\mathbf{n}} -\\left(\\frac{\\partial \\ln \\phi_i}{\\partial T}\\right)'_{p,\\mathbf{n}}\\right]$$\n",
      "    \n",
      "$$\\frac{\\partial F_i}{\\partial \\ln p} = p\\left[ \\left(\\frac{\\partial \\ln \\phi_i}{\\partial p}\\right)''_{T,\\mathbf{n}} -\\left(\\frac{\\partial \\ln \\phi_i}{\\partial p}\\right)'_{T,\\mathbf{n}}\\right]$$\n",
      "    \n",
      "$$ \\frac{\\partial F_i}{\\partial \\ln K_j} = \\frac{K_jz_j}{(1-\\beta+\\beta K_j)^2}[(1-\\beta)\\phi_{ij}''+\\beta\\phi_{ij}']+\\zeta $$\n",
      "\n",
      "$\\zeta$ is the Kronecker delta or $\\zeta = 0$ for $i\\neq j$ , and  $\\zeta = 0$ for $i=j$.  Also\n",
      "\n",
      "$$\\phi_{ij} = n\\left( \\frac{\\partial \\ln \\phi_i}{\\partial n_j}\\right)_{T,p}$$\n",
      "\n",
      "For the $F_{N+1}$ term,\n",
      "\n",
      "$$\\frac{\\partial F_{N+1}}{\\partial \\ln K_j}=\\frac{K_jz_j}{(1-\\beta+\\beta K_j)^2}$$\n",
      "\n",
      "and all other partials of $F_{N+1}$ in the Jacobian are zero.  For the specified term\n",
      "\n",
      "$$\\frac{\\partial F_{N+2}}{X_s}=1$$\n",
      "    \n",
      "and all other partials of $F_{N+2}$ in the Jacobian are zero.\n",
      "\n",
      "From GERG 2004 Monograph, Eqn 7.27:\n",
      "\n",
      "\n",
      "$$\\ln \\phi_i  = \\left( \\frac{\\partial n\\alpha^r}{\\partial n_i}\\right)_{T,V,n_j}-\\ln Z$$\n",
      "        \n",
      "and (Kunz, 2012, Table B4)\n",
      "\n",
      "$$\\left( \\frac{\\partial n\\alpha^r}{\\partial n_i}\\right)_{T,V,n_j} = \\alpha^r + n\\left( \\frac{\\partial \\alpha^r}{\\partial n_i}\\right)_{T,V,n_j}$$\n",
      "        \n",
      "so\n",
      "\n",
      "$$\\ln \\phi_i  = \\alpha^r + n\\left( \\frac{\\partial \\alpha^r}{\\partial n_i}\\right)_{T,V,n_j}-\\ln Z$$"
     ]
    },
    {
     "cell_type": "heading",
     "level": 2,
     "metadata": {},
     "source": [
      "Density marching phase envelope construction(T,P)"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Another two alternatives have been proposed in A DENSITY MARCHING METHOD FOR CALCULATING PHASE ENVELOPES\n",
      "Gadhiraju Venkatarathnam, I&ECR, 2014\n",
      "\n",
      "In this paper, density marching methods are proposed rather than methods that march in temperature, pressure, or K-factor.\n",
      "\n",
      "The system of equations to be solved is similar to that of the GERG 2004 formulation, where the unknowns are $\\ln(T)$, $\\ln(p)$, and $\\ln(K_i)$\n",
      "\n",
      "(A1) - OK\n",
      "$$F_i = \\ln K_i+\\ln\\phi(T,p,\\mathbf{y})-\\ln \\phi(T,p,\\mathbf{x})=0,  i=1,2,3... N$$\n",
      "\n",
      "(A2) - OK\n",
      "$$F_{N+1} = \\sum_{i=1}^{N}\\frac{z_i(K_i-1)}{1-\\beta+\\beta K_i}=0$$\n",
      "    \n",
      "(A3) - TYPO, should be $\\ln(\\rho)$ rather than $\\rho$, and should be all on left-hand-side\n",
      "$$F_{N+2} = \\ln \\rho ''-\\ln\\rho''_{S} = 0$$\n",
      "\n",
      "(A6) - TYPO, missing an n to multiply the terms $\\left( \\frac{\\partial \\ln \\phi_i}{\\partial n_j}\\right)_{T,p}$\n",
      "$$ \\frac{\\partial F_i}{\\partial \\ln K_j} = \\frac{K_jz_j}{(1-\\beta+\\beta K_j)^2}[(1-\\beta)\\phi_{ij}''+\\beta\\phi_{ij}']+\\zeta $$\n",
      "\n",
      "$\\zeta$ is the Kronecker delta or $\\zeta = 0$ for $i\\neq j$ , and  $\\zeta = 0$ for $i=j$.\n",
      "\n",
      "(A7) - OK\n",
      "$$\\frac{\\partial F_i}{\\partial \\ln T} = T\\left[ \\left(\\frac{\\partial \\ln \\phi_i}{\\partial T}\\right)''_{p,\\mathbf{n}} -\\left(\\frac{\\partial \\ln \\phi_i}{\\partial T}\\right)'_{p,\\mathbf{n}}\\right]$$\n",
      "\n",
      "(A8) - OK\n",
      "$$\\frac{\\partial F_i}{\\partial \\ln p} = p\\left[ \\left(\\frac{\\partial \\ln \\phi_i}{\\partial p}\\right)''_{T,\\mathbf{n}} -\\left(\\frac{\\partial \\ln \\phi_i}{\\partial p}\\right)'_{T,\\mathbf{n}}\\right]$$\n",
      "    \n",
      "(A9) - OK\n",
      "$$\\frac{\\partial F_{N+1}}{\\partial \\ln K_j}=\\frac{K_jz_j}{(1-\\beta+\\beta K_j)^2}$$\n",
      "\n",
      "(A11) - OK\n",
      "$$ \\frac{\\partial F_{N+2}}{\\partial \\ln K_j} = \\frac{K_jz_j(1-\\beta)\\beta}{(1-\\beta+\\beta K_j)^2}\\left(n\\left(\\frac{\\partial \\rho}{\\partial n_j}\\right)''_{T,p}\\right)\\left(\\frac{1}{\\rho''_{S}}\\right)$$\n",
      "\n",
      "(A12) - OK\n",
      "$$\\frac{\\partial F_{N+2}}{\\partial \\ln T}=\\left(\\frac{\\partial \\rho}{\\partial T}\\right)''_{p,n}\\frac{T}{\\rho''_{S}}$$\n",
      "\n",
      "(A13) - OK\n",
      "$$\\frac{\\partial F_{N+2}}{\\partial \\ln p}=\\left(\\frac{\\partial \\rho}{\\partial p}\\right)''_{T,n}\\frac{p}{\\rho''_{S}}$$\n"
     ]
    },
    {
     "cell_type": "heading",
     "level": 2,
     "metadata": {},
     "source": [
      "Density marching phase envelope construction (density)"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Another two alternatives have been proposed in A DENSITY MARCHING METHOD FOR CALCULATING PHASE ENVELOPES\n",
      "Gadhiraju Venkatarathnam, I&ECR, 2014\n",
      "\n",
      "In this paper, density marching methods are proposed rather than methods that march in temperature, pressure, or K-factor.\n",
      "\n",
      "The system of equations to be solved is similar to that of the GERG 2004 formulation, where the unknowns are $\\ln(T)$, $\\ln(p)$, and $\\ln(K_i)$\n",
      "\n",
      "(A14) - OK\n",
      "$$F_i = \\ln K_i+\\ln\\phi(T,\\rho'',\\mathbf{y})-\\ln \\phi(T,\\rho',\\mathbf{x})=0,  i=1,2,3... N$$\n",
      "\n",
      "(A15) - OK\n",
      "$$F_{N+1} = \\sum_{i=1}^{N}\\frac{z_i(K_i-1)}{1-\\beta+\\beta K_i}=0$$\n",
      "    \n",
      "(A16) - OK\n",
      "$$F_{N+2} = p(T,\\rho',\\mathbf{y})-p(T,\\rho'',\\mathbf{x}) = 0$$\n",
      "\n",
      "(A17) - TYPO, missing an n to multiply the terms $\\left( \\frac{\\partial \\ln \\phi_i}{\\partial n_j}\\right)_{T,p}$\n",
      "$$ \\frac{\\partial F_i}{\\partial \\ln K_j} = \\frac{K_jz_j}{(1-\\beta+\\beta K_j)^2}[(1-\\beta)\\phi_{ij}''+\\beta\\phi_{ij}']+\\zeta $$\n",
      "\n",
      "$\\zeta$ is the Kronecker delta or $\\zeta = 0$ for $i\\neq j$ , and  $\\zeta = 0$ for $i=j$.\n",
      "\n",
      "(A18) - OK\n",
      "$$\\frac{\\partial F_i}{\\partial \\ln T} = T\\left[ \\left(\\frac{\\partial \\ln \\phi_i}{\\partial T}\\right)''_{p,\\mathbf{n}} -\\left(\\frac{\\partial \\ln \\phi_i}{\\partial T}\\right)'_{p,\\mathbf{n}}\\right]$$\n",
      "\n",
      "(A19) - OK\n",
      "$$\\frac{\\partial F_i}{\\partial \\ln \\rho'} = -\\rho'\\left(\\frac{\\partial \\ln \\phi_i}{\\partial \\rho}\\right)'_{T,n}$$\n",
      "    \n",
      "(A20) - OK\n",
      "$$\\frac{\\partial F_{N+1}}{\\partial \\ln K_j}=\\frac{K_jz_j}{(1-\\beta+\\beta K_j)^2}$$\n",
      "\n",
      "(A22) - TYPO Second derivative of ln(phi) with respect to rho' needs constraints, first needs to have the constraints in the right place\n",
      "$$ \\frac{\\partial F_{N+2}}{\\partial \\ln K_j} = \\frac{RTK_jz_j}{(1-\\beta+\\beta K_j)^2}\\left[(1-\\beta)\\left(\\frac{\\partial \\ln \\phi_i}{\\partial \\rho}\\right)''_{T,n}+\\beta\\left(\\frac{\\partial \\ln \\phi_i}{\\partial \\rho}\\right)'_{T,n}\\right]$$\n",
      "\n",
      "(A23) - TYPO Should be A23\n",
      "$$\\frac{\\partial F_{N+2}}{\\partial \\ln T}=T\\left[\\left(\\frac{\\partial p}{\\partial T}\\right)'_{\\rho',n}-\\left(\\frac{\\partial p}{\\partial T}\\right)''_{\\rho'',n} \\right]$$\n",
      "\n",
      "(A24) - TYPO Should be A24\n",
      "$$\\frac{\\partial F_{N+2}}{\\partial \\ln p'}=\\rho'\\left(\\frac{\\partial p}{\\partial \\rho}\\right)'_{T,n}$$\n"
     ]
    },
    {
     "cell_type": "heading",
     "level": 2,
     "metadata": {},
     "source": [
      "Other analytic derivatives"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Three analyic derivatives are not provided in GERG and need to be rederived:\n",
      "\n",
      "$$\\left(\\frac{\\partial \\ln \\phi_i}{\\partial \\rho}\\right)_{T,n}$$\n",
      "\n",
      "$$\\left(\\frac{\\partial \\ln \\phi_i}{\\partial T}\\right)_{\\rho,n}$$\n",
      "\n",
      "$$n\\left(\\frac{\\partial \\rho}{\\partial n_j}\\right)_{T,p}$$\n",
      "The last is for T,p marching, the first two are for density marching."
     ]
    },
    {
     "cell_type": "heading",
     "level": 3,
     "metadata": {},
     "source": [
      "Derivations"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "For $n\\left( \\frac{\\partial \\rho}{\\partial n_i}\\right)_{T,p,n_j}$\n",
      "-----------------------------------------------------------------\n",
      "\n",
      "GERG 2007 Monograph Equation 7.32 gives\n",
      "\n",
      "$$\\left( \\frac{\\partial V}{\\partial n_i}\\right)_{T,p,n_j} = \\dfrac{-\\left(\\dfrac{\\partial p}{\\partial n_i}\\right)_{T,V,n_j}}{\\left(\\dfrac{\\partial p}{\\partial V}\\right)_{T,n}}$$\n",
      "\n",
      "expand the left hand side with \n",
      "\n",
      "$V = vn = \\dfrac{n}{\\rho}$\n",
      "\n",
      "n held constant in derivative, so get\n",
      "\n",
      "$$\\left( \\frac{\\partial V}{\\partial n_i}\\right)_{T,p,n_j} = n\\left( \\frac{\\partial (1/\\rho)}{\\partial n_i}\\right)_{T,p,n_j} = -\\frac{n}{\\rho^2}\\left( \\frac{\\partial \\rho}{\\partial n_i}\\right)_{T,p,n_j}$$\n",
      "\n",
      "so\n",
      "\n",
      "$$ n\\left( \\frac{\\partial \\rho}{\\partial n_i}\\right)_{T,p,n_j} = -\\rho^2\\left( \\frac{\\partial V}{\\partial n_i}\\right)_{T,p,n_j} $$\n",
      "\n",
      "For $\\left(\\frac{\\partial \\ln \\phi_i}{\\partial \\rho}\\right)_{T,n}$ and $\\left(\\frac{\\partial \\ln \\phi_i}{\\partial T}\\right)_{\\rho,n}$\n",
      "--------------------------------------------------------------------------------------------------------------------------------------------\n",
      "\n",
      "GERG 2007 Monograph 7.27\n",
      "$$\\ln \\phi_i  = \\left(\\frac{\\partial n \\alpha^r}{\\partial n_i} \\right)_{T,V,n_j} - \\ln Z $$\n",
      "\n",
      "GERG 2007 Monograph 7.34\n",
      "$$\\ln \\left(\\frac{f_i}{n_i}\\right) = \\ln\\left(\\frac{RT}{V}\\right)+\\left(\\frac{\\partial n\\alpha^r}{\\partial n_i} \\right)_{T,V,n_j}$$\n",
      "\n",
      "and $Z = (pV)/(nRT)$, thus\n",
      "\n",
      "$$\\ln \\phi_i = \\ln \\left(\\frac{f_i}{n_i}\\right) - \\ln\\left(\\frac{RT}{V}\\right)  - \\ln \\left(\\frac{pV}{nRT}\\right) ?=  \\ln \\left(\\frac{f_i}{n_i}\\right) -\\ln\\left(\\frac{p}{n}\\right)$$\n",
      "\n",
      "$$\\frac{\\partial \\left[\\ln\\left(\\frac{p}{n}\\right)\\right]}{\\partial T} = \\frac{1}{pn}\\left(\\frac{\\partial p}{\\partial T}\\right)_{\\rho,n}$$\n"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Equivalent fugacities for the $i$-th component\n",
      "\n",
      "$$ F_k^A = \\ln f_i(T,p,\\mathbf{x})-\\ln f_i(T,p,\\mathbf{y}) = 0\\mbox{ for } k = i = 1...N $$\n",
      "\n",
      "Material balance\n",
      "\n",
      "$$ F_k^B = \\frac{z_i-x_i}{y_i-x_i}-\\frac{z_{N-1}-x_{N-1}}{y_{N-1}-x_{N-1}}\\mbox{  for  }i=1..N-2; k = i+N; k = N+1..2N-2$$\n",
      "\n",
      "The independent variables to be obtained are\n",
      "\n",
      "$$ \n"
     ]
    },
    {
     "cell_type": "heading",
     "level": 1,
     "metadata": {},
     "source": [
      "Conversion of derivatives"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "To convert partial with $T$, $V$, $x_k$ held constant to one with $\\tau$, $\\delta$, $x_k$ held constant, use Gernert 3.118, or\n",
      "\n",
      "$$ \\frac{\\partial}{\\partial x_j} [Y]_{T,V,x_k} = \\frac{\\partial}{\\partial x_j} [Y]_{\\tau,\\delta,x_k}+\\left(\\frac{\\partial\\delta}{\\partial x_j}\\right)_{T,V,x_k}\\left.\\frac{\\partial Y}{\\partial\\delta}\\right|_{\\tau,\\bar x}+\\left(\\frac{\\partial\\tau}{\\partial x_j}\\right)_{T,V,x_k}\\left.\\frac{\\partial Y}{\\partial\\tau}\\right|_{\\delta,\\bar x} $$\n",
      "\n",
      "To convert pressure,\n",
      "\n",
      "$$ p=\\rho R T(1+\\delta \\alpha_\\delta) $$\n",
      "\n",
      "convert $\\rho$ and $T$ to reduced variables\n",
      "\n",
      "$$ p=\\rho_r(\\bar x)\\delta R \\frac{T_r(\\bar x)}{\\tau}(1+\\delta \\alpha_\\delta) = \\rho_r(\\bar x)R \\frac{T_r(\\bar x)}{\\tau}\\delta (1+\\delta \\alpha_\\delta)$$\n",
      "\n",
      "All the derivatives\n",
      "\n",
      "$$ \\frac{dp}{d\\tau}\\times\\frac{d\\tau}{dx_j} = -\\rho_r(\\bar x)\\delta R \\frac{T_r(\\bar x)}{\\tau^2}(1+\\delta \\alpha_\\delta) \\times \\frac{1}{T}\\frac{\\partial T_r}{\\partial x_j}|_{T,V,x_k}$$\n",
      "\n",
      "$$ \\frac{dp}{d\\delta}\\times\\frac{d\\delta}{dx_j} = \\rho_r(\\bar x) R \\frac{T_r(\\bar x)}{\\tau}[ (1+\\delta \\alpha_\\delta)+ \\delta(\\alpha_\\delta+\\delta \\alpha_{\\delta\\delta})] \\times \\frac{-\\delta}{\\rho_r}\\frac{\\partial \\rho_r}{\\partial x_j}|_{T,V,x_k}$$\n",
      "\n",
      "$$ \\frac{\\partial p}{\\partial x_j}|_{\\tau,\\delta,x_k} = \\frac{\\delta R}{\\tau}\\left[\\rho_r(\\bar x)T_r(\\bar x)(\\delta \\frac{\\partial}{\\partial x_j}[\\alpha_\\delta]_{\\tau, \\delta, x_k})+(1+\\delta \\alpha_\\delta)\\left(\\rho_r \\frac{\\partial T_r}{\\partial x_j}+T_r \\frac{\\partial \\rho_r}{\\partial x_j}\\right)\\right]$$\n",
      "\n",
      "with $\\delta$ factored out"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Turn back into normal variables\n",
      "$$ \\frac{dp}{d\\tau}\\times\\frac{d\\tau}{dx_j} = -\\frac{\\rho R}{\\tau}(1+\\delta \\alpha_\\delta)\\frac{\\partial T_r}{\\partial x_j}|_{T,V,x_k}$$\n",
      "\n",
      "$$ \\frac{dp}{d\\delta}\\times\\frac{d\\delta}{dx_j} = - R T[ (1+\\delta \\alpha_\\delta)+ \\delta(\\alpha_\\delta+\\delta \\alpha_{\\delta\\delta})]\\delta\\frac{\\partial \\rho_r}{\\partial x_j}|_{T,V,x_k}$$\n",
      "\n",
      "$$ \\frac{\\partial p}{\\partial x_j}|_{\\tau,\\delta,x_k} = R\\left[\\rho T (\\delta \\frac{\\partial}{\\partial x_j}[\\alpha_\\delta]_{\\tau, \\delta, x_k})+(1+\\delta \\alpha_\\delta)\\left(\\frac{\\rho}{\\tau} \\frac{\\partial T_r}{\\partial x_j}+\\delta T \\frac{\\partial \\rho_r}{\\partial x_j}\\right)\\right]$$"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "First term cancels with a term in the third one, yielding\n",
      "\n",
      "$$ \\frac{dp}{d\\delta}\\times\\frac{d\\delta}{dx_j} = -\\delta R T[ (1+\\delta \\alpha_\\delta)+ \\delta(\\alpha_\\delta+\\delta \\alpha_{\\delta\\delta})]\\frac{\\partial \\rho_r}{\\partial x_j}|_{T,V,x_k}$$\n",
      "\n",
      "$$ \\frac{\\partial p}{\\partial x_j}|_{\\tau,\\delta,x_k} = \\rho R T (\\delta \\frac{\\partial}{\\partial x_j}[\\alpha_\\delta]_{\\tau, \\delta, x_k})+(1+\\delta \\alpha_\\delta)\\delta R T \\frac{\\partial \\rho_r}{\\partial x_j}$$"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "First term in first line cancels with term at end of second line, yielding\n",
      "\n",
      "$$ \\frac{dp}{d\\delta}\\times\\frac{d\\delta}{dx_j} = -\\delta R T[ \\delta(\\alpha_\\delta+\\delta \\alpha_{\\delta\\delta})]\\frac{\\partial \\rho_r}{\\partial x_j}|_{T,V,x_k}$$\n",
      "\n",
      "$$ \\frac{\\partial p}{\\partial x_j}|_{\\tau,\\delta,x_k} = \\rho R T (\\delta \\frac{\\partial}{\\partial x_j}[\\alpha_\\delta]_{\\tau, \\delta, x_k})$$"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Total equation is\n",
      "\n",
      "$$ \\frac{\\partial p}{\\partial x_j}|_{T,V,x_k} = \\rho R T \\left(\\delta \\frac{\\partial}{\\partial x_j}[\\alpha_\\delta]_{\\tau, \\delta, x_k}-\\frac{\\delta}{\\rho_r}(\\alpha_\\delta+\\delta \\alpha_{\\delta\\delta})\\frac{\\partial \\rho_r}{\\partial x_j}|_{T,V,x_k}\\right)$$"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "$$ p=\\rho R T(1+\\delta \\alpha_\\delta) $$\n",
      "\n",
      "$$ p=\\rho R T(1+\\frac{\\rho}{\\rho_r} \\alpha_\\delta) $$"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "$$ \\frac{\\partial p}{\\partial x_j}|_{T,V,x_k} = \\rho R T \\left(\\rho\\frac{-1}{\\rho_r^2}\\left(\\frac{\\partial \\rho_r}{\\partial x_j}\\right)_{x_{k\\neq j}}\\alpha_\\delta + \\frac{\\rho}{\\rho_r}\\left(\\frac{\\partial}{\\partial x_j}\\left(\\frac{\\partial \\alpha_r}{\\partial \\delta}_{\\tau,\\bar x}\\right)\\right)_{T,V,x_{k\\neq j}}\\right)$$\n",
      "\n",
      "$$ \\frac{\\partial p}{\\partial x_j}|_{T,V,x_k} = \\delta\\rho R T \\left(\\frac{-1}{\\rho_r}\\left(\\frac{\\partial \\rho_r}{\\partial x_j}\\right)_{x_{k\\neq j}}\\alpha_\\delta + \\left(\\frac{\\partial}{\\partial x_j}\\left(\\frac{\\partial \\alpha_r}{\\partial \\delta}_{\\tau,\\bar x}\\right)\\right)_{T,V,x_{k\\neq j}}\\right)$$"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import numpy as np\n",
      "import matplotlib.pyplot as plt\n",
      "import CoolProp.CoolProp as CP\n",
      "%matplotlib inline"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 1
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "p = np.linspace(1000, 12000, 100)\n",
      "mix = 'HEOS::Water[0.7]&Ethanol[0.3]'\n",
      "rhoL = CP.PropsSI('D','P',p,'Q',[0]*len(p),mix)\n",
      "rhoV = CP.PropsSI('D','P',p,'Q',[1]*len(p),mix)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 6
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plt.plot(rhoL,p,rhoV,p)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 7,
       "text": [
        "[<matplotlib.lines.Line2D at 0xa9fe908>,\n",
        " <matplotlib.lines.Line2D at 0xa9feb38>]"
       ]
      },
      {
       "metadata": {},
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEACAYAAACpoOGTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAD5RJREFUeJzt3G+MVNd5x/HvlMVKIiumpBUYlgoXYxVSNSVtyEatyjiu\nJUJTcBWpFlJFYks1aktitakDxFK9fRXbqWrLQrGRi6JNmoi61LJIg/9Qy6O+aYid2OvGsDaLQwpE\nXkeJ/CKWqmAxfXEO3unwLDu75w4zDN+PNNr759w7z9Gw97f3nDuAJEmSJEmSJEmSJEmSJEmSNJA2\nAhPAcWDnDG0eyvvHgXUt2xcBB4BjwFFgpHtlSpIupQXAJLASWAi8BKxpa7MJOJSXPwp8p2XfGHB7\nXh4CrulWoZKkS+tjwFMt67vyq9UjwK0t6xPAElIYvN7V6iRJ8/ZLhccvB061rJ/O22ZrMwxcB/wE\n+CrwfeBR4H2F9UiSKlIaEM0O29WC44aADwNfyT/f5sK7D0lSjwwVHn8GWNGyvoJ0h3CxNsN5Wy23\nfT5vP0AQEKtWrWqeOHGisExJuuKcAK4vOUHpHcQLwGrSJPVVpLmGg21tDgLb8vII8BYwBbxBGnq6\nIe/7Q+CV9jc4ceIEzWZzYF/33HNPz2uwf/bvSuvbldA/YFXh9b34DuIdYAfwNOmJpn2kR1a35/17\nSU8wbSI97fQ2cFvL8Z8FvkEKlxNt+yRJPVQaEABP5lervW3rO2Y4dhz4SAU1SJIqVjrEpEL1er3X\nJXSV/bt8DXLfYPD7V4X2p4v6UTOPp0mSOlSr1aDwGu8dhCQpZEBIkkIGhCQpZEBIkkIGhCQpZEBI\nkkIGhCQpZEBIkkIGhCQpZEBIkkIGhCQpZEBIkkIGhCQpZEBIkkIGhCQpZEBIkkIGhCQpZEBIkkIG\nhCQpZEBIkkIGhCQpZEBIkkIGhCQpZEBIkkIGhCQpZEBIkkJVBMRGYAI4Duycoc1Def84sK5t3wLg\nReBbFdQiSapIaUAsAPaQQmItsBVY09ZmE3A9sBq4A3i4bf+dwFGgWViLJKlCpQGxHpgETgJngf3A\nlrY2m4GxvHwEWAQsyevDpAD5J6BWWIskqUKlAbEcONWyfjpv67TNA8BdwLnCOiRJFSsNiE6Hhdrv\nDmrAJ4E3SfMP3j1IUp8ZKjz+DLCiZX0F6Q7hYm2G87ZPkYafNgHvAd4PfA3Y1v4mo6Oj7y7X63Xq\n9Xph2ZI0WBqNBo1Go9Jzlv7lPgS8CtwE/Bj4Lmmi+lhLm03AjvxzBHgw/2y1Afhb4I+D92g2m85f\nS9Jc1Go1KLzGl95BvEO6+D9NeqJpHykctuf9e4FDpHCYBN4GbpvhXKaAJPWRy2Hs3zsISZqjKu4g\n/Ca1JClkQEiSQgaEJClkQEiSQgaEJClkQEiSQgaEJClkQEiSQgaEJClkQEiSQgaEJClkQEiSQgaE\nJClkQEiSQgaEJClkQEiSQgaEJClkQEiSQgaEJClkQEiSQgaEJClkQEiSQgaEJClkQEiSQgaEJClk\nQEiSQgaEJClkQEiSQgaEJClURUBsBCaA48DOGdo8lPePA+vythXAc8ArwA+Az1VQiySpIqUBsQDY\nQwqJtcBWYE1bm03A9cBq4A7g4bz9LPDXwAeBEeCvgmMlST1SGhDrgUngJOmCvx/Y0tZmMzCWl48A\ni4AlwBvAS3n7z4FjwLLCeiRJFSkNiOXAqZb103nbbG2G29qsJA09HSmsR5JUkaHC45sdtqtd5Lir\ngQPAnaQ7iQuMjo6+u1yv16nX6x0XKElXgkajQaPRqPSc7RfuuRoBRklzEAC7gXPAfS1tHgEapOEn\nSBPaG4ApYCHw78CTwIMzvEez2ew0hyRJALVaDQqv8aVDTC+QJp9XAlcBtwIH29ocBLbl5RHgLVI4\n1IB9wFFmDgdJUo+UDjG9A+wAniY90bSPNNm8Pe/fCxwiPck0CbwN3Jb3/R7wZ8DLwIt5227gqcKa\nJEkVKB1iuhQcYpKkOeqHISZJ0oAyICRJIQNCkhQyICRJIQNCkhQyICRJIQNCkhQyICRJIQNCkhQy\nICRJIQNCkhQyICRJIQNCkhQyICRJIQNCkhQyICRJIQNCkhQyICRJIQNCkhQyICRJIQNCkhQyICRJ\nIQNCkhQyICRJIQNCkhQyICRJIQNCkhQyICRJoSoCYiMwARwHds7Q5qG8fxxYN8djJUk9UBoQC4A9\npAv9WmArsKatzSbgemA1cAfw8ByOlST1SGlArAcmgZPAWWA/sKWtzWZgLC8fARYBSzs8VpLUI6UB\nsRw41bJ+Om/rpM2yDo6VJPXIUOHxzQ7b1UreZHR09N3ler1OvV4vOZ0kDZxGo0Gj0aj0nEUXbmAE\nGCXNIwDsBs4B97W0eQRokIaQIE1KbwCu6+BYgGaz2WkOSZIAarUaFF7jS4eYXiBNPq8ErgJuBQ62\ntTkIbMvLI8BbwFSHx0qSeqR0iOkdYAfwNOmppH3AMWB73r8XOER6kmkSeBu4bZZjJUl9oHSI6VJw\niEmS5qgfhpgkSQPKgJAkhQwISVLIgJAkhQwISVLIgJAkhQwISVLIgJAkhQwISVLIgJAkhQwISVLI\ngJAkhQwISVLIgJAkhQwISVLIgJAkhQwISVLIgJAkhQwISVLIgJAkhQwISVLIgJAkhQwISVLIgJAk\nhQwISVLIgJAkhQwISVLIgJAkhUoDYjFwGHgNeAZYNEO7jcAEcBzY2bL9y8AxYBx4HLimsB5JUkVK\nA2IXKSBuAJ7N6+0WAHtIIbEW2AqsyfueAT4IfIgUMrsL65EkVaQ0IDYDY3l5DLglaLMemAROAmeB\n/cCWvO8wcC4vHwGGC+uRJFWkNCCWAFN5eSqvt1sOnGpZP523tbsdOFRYjySpIkMdtDkMLA223922\n3syvdtG26Fy/AL4Z7RwdHX13uV6vU6/XOzilJF05Go0GjUaj0nPWCo+fAOrAG8C1wHPAb7S1GQFG\nSXMQkOYZzgH35fXPAH8O3AT8b/AezWazk4yRJJ1Xq9Wg8BpfOsR0EPh0Xv408ETQ5gVgNbASuAq4\nNR8HKTTuIs1JROEgSeqR0juIxcBjwK+RJqH/FHgLWAY8CvxRbvcJ4EHSE037gC/l7cdJofGzvP5f\nwF+2vYd3EJI0R1XcQZQGxKVgQEjSHPXDEJMkaUAZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZ\nEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKk\nkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkUElALAYOA68BzwCLZmi3EZgAjgM7\ng/2fB87l80mS+kRJQOwiBcQNwLN5vd0CYA8pJNYCW4E1LftXADcDPyqoQ5LUBSUBsRkYy8tjwC1B\nm/XAJHASOAvsB7a07P9H4AsFNUiSuqQkIJYAU3l5Kq+3Ww6calk/nbdBCorTwMsFNUiSumRolv2H\ngaXB9rvb1pv51S7aBvBe4Iuk4aXzarPUIkm6hGYLiJsvsm+KFB5vANcCbwZtzpDmGc5bQbprWAWs\nBMbz9mHge6QhqQvOMzo6+u5yvV6nXq/PUrYkXVkajQaNRqPSc5b81X4/8FPgPtIE9SIunKgeAl4F\nbgJ+DHyXNFF9rK3dD4HfAX4WvE+z2ZzpRkSSFKnValA4MlMyB3Ev6Q7jNeDjeR1gGfDtvPwOsAN4\nGjgK/AsXhgPMPBQlSeqRy2Hc3zsISZqjXt9BSJIGmAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKk\nkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEh\nSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQoZEJKkkAEhSQqVBMRi4DDwGvAMsGiGdhuBCeA4\nsLNt32eBY8APgPsKapEkVawkIHaRAuIG4Nm83m4BsIcUEmuBrcCavO9GYDPwW8BvAv9QUMtlq9Fo\n9LqErrJ/l69B7hsMfv+qUBIQm4GxvDwG3BK0WQ9MAieBs8B+YEve9xfAl/J2gJ8U1HLZGvR/pPbv\n8jXIfYPB718VSgJiCTCVl6fyervlwKmW9dN5G8Bq4A+A7wAN4HcLapEkVWxolv2HgaXB9rvb1pv5\n1S7a1vrevwyMAB8BHgN+fZZ6JEmXgQmmw+PavN5uBHiqZX030xPVTwIbWvZNAh8IzjHJdAD58uXL\nl6/OXpP00P1MX+x3AfcGbYaAE8BK4CrgJaYnqbcDf5+XbwD+p1uFSpIurcXAf3DhY67LgG+3tPsE\n8CopzXa3bF8IfB34b+B7QL275UqSJEkaGIP+pbsq+gfweeBcPl8/Ke3fl0mf3TjwOHBN1yrt3Gyf\nBcBDef84sG6Ox/bafPu3AngOeIX0u/a57pY5byWfH6TvcL0IfKtbBRYo6dsi4ADp9+0oaZ64790P\nfCEv7ySez1hAGqZaSRqeap3PuJF0gVqY13+1W4XOU2n/IP1iPgX8kP4LiNL+3cz0I9f3znD8pTTb\nZwGwCTiUlz9Kely702N7raR/S4HfzstXk4aPB6l/5/0N8A3gYNeqnJ/Svo0Bt+flIfrjj7FZTTD9\nPYqlxE9EfYz//0TULqa/vf0Y8PGuVVeutH8A/0r61nk/BkQV/TvvT4B/rrS6ueuk1keAW1vWzz/V\n12k/e2m+/Yu+6/QEcFOl1ZUr7d8waX71RvrvDqKkb9cAr8/lzfrlP+sb9C/dlfZvS15/uVsFFirt\nX6vbmf7rp1c6qXWmNss6OLbX5tu/4bY2K0nDF0cqrq9UyecH8ABwF2k4t9+UfHbXkf7Hiq8C3wce\nBd53sTeb7YtyVRr0L911q3/vBb5IGoY5rzbn6sp18/NrPdcvgG/OrbTKdVIr9OZzqMJ8+9d63NWk\nsew7gZ9XUVSF5tu/GvBJ4E3S/EO9wpqqUvLZDQEfBnYAzwMPku4+/m6mk1zKgLj5IvumSBefN0hf\nunszaHOGNA5/3gpSMpJ/Pp6Xnycl/weAnxbUO1fd6t8q0l9q43n7MOmx4PUznKdbuvn5AXyGNHba\nD8MVs9UatRnObRZ2cGyvzbd/Z/LyQuDfSEOBT3SpxhIl/fsU6f+Z2wS8B3g/8DVgW7eKnaOSvtVy\n2+fz9gP03/BnaNC/dFfav1b9OAdR2r+NpKdifqWrVXauk8+idSJwhOmJwE4/x14q6V+NdMF8oOtV\nzl9J/1ptoP/mIEr79p+kayTAKP33xGdo0L90V9q/Vq/TfwFR2r/jwI9It/UvAl/pcr2diGrdnl/n\n7cn7x0m37hc7tt/Mt3+/T7pDf4npz2vjJah3rko+v/M20H9PMUFZ3z5EuoPop0fKJUmSJEmSJEmS\nJEmSJEmSJEmSJEnqjv8DwHN/hNOHRlAAAAAASUVORK5CYII=\n",
       "text": [
        "<matplotlib.figure.Figure at 0xaa7e550>"
       ]
      }
     ],
     "prompt_number": 7
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "$$\\ln f_i = \\ln (x_i p \\phi_i)$$\n",
      "\n",
      "$$\\frac{\\partial \\ln f_i}{\\partial x_j} = \\frac{\\partial \\ln (x_i p \\phi_i)}{\\partial x_j} = \\frac{p}{x_i p \\phi_i}\n",
      "\\left(\\frac{\\partial x_i}{\\partial x_j}\\phi_i + x_i \\frac{\\partial x_i}{\\partial x_j}\\right) = \\frac{1}{x_i \\phi_i}\n",
      "\\left(\\frac{\\partial x_i}{\\partial x_j}\\phi_i + x_i \\frac{\\partial x_i}{\\partial x_j}\\right) = \n",
      "\\frac{\\partial x_i}{\\partial x_j}\\frac{1}{x_i} + \\frac{1}{\\phi_i} \\frac{\\partial \\phi_i}{\\partial x_j}$$\n",
      "\n",
      "$$\\frac{\\partial \\ln f_i}{\\partial x_j} = \\frac{\\partial \\ln (x_i p \\phi_i)}{\\partial x_j} = \n",
      "\\frac{1}{x_i}\\frac{\\partial x_i}{\\partial x_j} + \\frac{\\partial \\ln\\phi_i}{\\partial x_j}$$\n",
      "\n",
      "Thus if $i=j$, $\\displaystyle \\frac{\\partial x_i}{\\partial x_j} = 1$, otherwise first term goes away\n",
      "\n",
      "If $i = N-1$, $\\displaystyle\\frac{\\partial x_{N-1}}{\\partial x_j} = -1$, since \n",
      "$$x_{N-1} = 1-\\sum_{k=1}^{N-1}x_k$$"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}