Merge pull request #129 from ipqa-research/uniquac_more_test

Uniquac more test

doc/page/usage/excessmodels/unifaclv.md

@@ -2,6 +2,8 @@
 title: UNIFAC-LV
 ---
 
+[TOC]
+
 # UNIFAC
 
 [[UNIFAC]] (UNIQUAC Functional-group Activity Coefficients) is an Excess Gibbs
@@ -15,25 +17,229 @@ $$ \frac{G^E}{RT} = \frac{G^{E,r}}{RT} + \frac{G^{E,c}}{RT} $$
 
 Being:
 
-- Combinatorial: Accounts for the size and shape of the molecules. The involved
-equations can be checked in the API documentation:
-[[yaeos__models_ge_group_contribution_unifac:Ge_combinatorial]]
+- Combinatorial: Accounts for the size and shape of the molecules.
 
-- Residual: Accounts for the energy interactions between different functional
-groups. The involved equations can be checked in the API documentation:
-[[yaeos__models_ge_group_contribution_unifac:Ge_residual]]
+- Residual: Accounts for the energy interactions between different functional groups.
 
 Each substance of a mixture modeled with [[UNIFAC]] must be represented by a
 list a functional groups and other list with the ocurrence of each functional
 group on the substance. The list of the functional groups culd be accesed on
 the DDBST web page:
 [https://www.ddbst.com/published-parameters-unifac.html](https://www.ddbst.com/published-parameters-unifac.html)
 
-# Examples
-## Calculating activity coefficients
-We can instantiate a [[UNIFAC]] model with a mixture  ethanol-water and
-evaluate the logarithm of activity coefficients of the model for a 0.5 mole
-fraction of each, and a temperature of 298.0 K.
+## Combinatorial term derivatives
+
+Calculate the UNIFAC combinatorial term of reduced Gibbs excess energy.
+The subroutine uses the Flory-Huggins and Staverman-Guggenheim
+combinatory terms as follows:
+
+### Flory-Huggins
+
+$$
+   G^{E,FH} =
+   RT \left(\sum_i^{NC} n_i \, \text{ln} \, r_i
+   - n \, \text{ln} \, \sum_j^{NC} n_j r_j
+   + n \, \text{ln} \, n \right)
+$$
+
+$$
+   \frac{dG^{E,FH}}{dn_i} =
+   RT \left(\text{ln} \, r_i - \text{ln} \, \sum_j^{NC} n_j r_j
+   + \text{ln} \, n + 1 - \frac{n r_i}{\displaystyle
+   \sum_j^{NC} n_j r_j} \right)
+$$
+
+$$
+   \frac{d^2G^{E,FH}}{dn_i dn_j} =
+   RT \left(- \frac{r_i + r_j}{\displaystyle \sum_l^{NC} n_l r_l}
+   + \frac{1}{n} + \frac{n r_i r_j}{\displaystyle \left(\sum_l^{NC}
+   n_l r_l \right)^2} \right)
+$$
+
+### Staverman-Guggenheim
+
+$$
+   \frac{G^{E,SG}}{RT} =
+   \frac{z}{2} \sum_i^{NC} n_i q_i
+   \left(\text{ln} \frac{q_i}{r_i}
+   - \text{ln} \, \sum_j^{NC} n_j q_j
+   + \text{ln} \, \sum_j^{NC} n_j r_j \right)
+$$
+
+$$
+   \frac{1}{RT}\frac{dG^{E,SG}}{dn_i} =
+   \frac{z}{2} q_i \left(
+   - \text{ln} \, \left(
+   \frac{r_i \sum_j^{NC} n_j q_j}{\displaystyle q_i \sum_j^{NC}
+   n_j r_j} \right) - 1 + \frac{\displaystyle r_i \sum_j^{NC} n_j
+   q_j}{\displaystyle q_i \sum_j^{NC} n_j r_j} \right)
+$$
+
+$$
+   \frac{1}{RT}\frac{d^2G^{E,SG}}{dn_i dn_j} =
+   \frac{z}{2} \left(- \frac{q_i q_j}{\displaystyle \sum_l^{NC} n_lq_l}
+   + \frac{q_i r_j + q_j r_i}{\displaystyle \sum_l^{NC} n_l r_l}
+   - \frac{\displaystyle r_i r_j \sum_l^{NC} n_l q_l}
+   {\left(\displaystyle \sum_l^{NC} n_l r_l \right)^2} \right)
+$$
+
+### Fredenslund et al. (UNIFAC)
+$$
+   \frac{G^{E,\text{UNIFAC}}}{RT} =
+   \frac{G^{E,FH}}{RT} + \frac{G^{E,SG}}{RT}
+$$
+
+## Residual terms derivatives
+Evaluate the UNIFAC residual term. The residual Gibbs excess energy
+and its derivatives are evaluated as:
+
+$$
+ \frac{G^{E,R}}{RT} = - \sum_i^{NC} n_i \sum_k^{NG} v_k^i Q_k
+ (\Lambda_k - \Lambda_k^i)
+$$
+
+With:
+
+$$
+ \Lambda_k = \text{ln} \, \sum_{j}^{NG} \Theta_j E_{jk}
+$$
+
+$$
+ \Lambda_k^i = \text{ln} \, \sum_{j}^{NG} \Theta_j^i E_{jk}
+$$
+
+$$
+ E_{jk} = \text{exp} \left(- \frac{U_{jk}}{RT} \right)
+$$
+
+$$
+ \Theta_j = \frac{Q_j \displaystyle \sum_{l}^{NC} n_l v_j^l}
+ {\displaystyle \sum_{k}^{NC} n_k \sum_{m}^{NG} v_m^l Q_m}
+$$
+
+$$
+ \Theta_j^i = \frac{Q_j v_j^i}{\displaystyle \sum_k^{NG} v_k^i Q_k}
+$$
+
+In the UNIFAC model, the \(\Theta_j^i \) values are calculated assuming
+that the molecule "i" is pure, hence only the subgroups of the molecule
+"i" must be considered for the calculation. On the other hand, for the
+\(\Theta_j \) values, all the system's subgroups are considered.
+
+### The compositional derivatives:
+
+$$
+ \frac{1}{R T} \frac{\partial G^{E,R}}{\partial n_\alpha} =
+ - \sum_k^{\mathrm{NG}} v_k^\alpha Q_k \left(\Lambda_k -
+ \Lambda_k^\alpha \right) - \sum_i^{\mathrm{NC}} n_i
+ \sum_k^{\mathrm{NG}} v_k^i Q_k
+ \frac{\partial \Lambda_k}{\partial n_\alpha}
+$$
+
+$$
+ \frac{1}{R T} \frac{\partial^2 G^{E,R}}{\partial n_
+ \alpha \partial n_\beta} = -\sum_k^{\mathrm{NG}} Q_k \left(v_k^\alpha
+ \frac{\partial \Lambda_k}{\partial n_\beta} + v_k^\beta
+ \frac{\partial \Lambda_k}{\partial n_\alpha}\right)
+ - \sum_k^{\mathrm{NG}} \left(\sum_i^{\mathrm{NC}} n_i v_k^i\right) Q_k
+ \frac{\partial^2 \Lambda_k}{\partial n_\alpha \partial n_\beta}
+$$
+
+With:
+
+$$
+ \frac{\partial \Lambda_k}{\partial n_\alpha}
+ = \frac{\sum_j^{\mathrm{NG}} v_j^\alpha Q_j E_{j k}}
+ {\sum_l^{\mathrm{NC}} n_l \sum_j^{\mathrm{NG}} v_j^l Q_j
+ E_{j k}} - \frac{\sum_m^{\mathrm{NG}} v_m^\alpha Q_m}
+ {\sum_l^{\mathrm{NC}} n_l \sum_m^{\mathrm{NG}} v_m^l Q_m}
+$$
+
+$$
+ \frac{\partial^2 \Lambda_k}{\partial n_\alpha \partial n_\beta}
+ = - \frac{\left(\sum_j^{\mathrm{NG}} v_j^\alpha Q_j E_{j k}\right)
+ \left(\sum_j^{\mathrm{NG}} v_j^\beta Q_j E_{j k}\right)}
+ {\left(\sum_l^{\mathrm{NC}} n_l \sum_j^{\mathrm{NG}} v_j^l Q_j
+ E_{j k}\right)^2} + \frac{\left(\sum_m^{\mathrm{NG}} v_m^\alpha
+ Q_m\right)\left(\sum_m^{\mathrm{NG}} v_m^\beta Q_m\right)}
+ {\left(\sum_l^{\mathrm{NC}} n_l
+ \sum_m^{\mathrm{NG}} v_m^l Q_m\right)^2}
+$$
+
+### The temperature derivatives:
+
+$$
+ \frac{\partial\left(\frac{G^{E, R}}{R T}\right)}{\partial T} =
+ -\sum_i^{\mathrm{NC}} n_i \sum_k^{\mathrm{NG}} v_k^i Q_k
+ \left(\frac{\partial \Lambda_k}{\partial T}
+ -\frac{\partial \Lambda_k^i}{\partial T}\right)
+$$
+
+$$
+ \frac{\partial^2\left(\frac{G^{E,R}}{R T}\right)}{\partial T^2} =
+ -\sum_i^{\mathrm{NC}} n_i \sum_k^{\mathrm{NG}} v_k^i Q_k
+ \left(\frac{\partial^2 \Lambda_k}{\partial T^2} -
+ \frac{\partial^2 \Lambda_k^i}{\partial T^2}\right)
+$$
+
+With:
+
+$$
+ \frac{\partial \Lambda_k}{\partial T} =
+ \frac{\sum_{j}^{NG} \Theta_j \frac{d E_{jk}}{dT}}
+ {\sum_{j}^{NG} \Theta_j E_{jk}}
+$$
+
+$$
+ \frac{\partial \Lambda_k^i}{\partial T} =
+ \frac{\sum_{j}^{NG} \Theta_j^i \frac{d E_{jk}}{dT}}
+ {\sum_{j}^{NG} \Theta_j^i E_{jk}}
+$$
+
+$$
+ \frac{\partial^2 \Lambda_k}{\partial T^2} =
+ \frac{\sum_{j}^{NG} \Theta_j \frac{d^2 E_{jk}}{dT^2}}
+ {\sum_{j}^{NG} \Theta_j E_{jk}}
+ - \left(\frac{\partial \Lambda_k}{\partial T} \right)^2
+$$
+
+$$
+ \frac{\partial^2 \Lambda_k^i}{\partial T^2} =
+ \frac{\sum_{j}^{NG} \Theta_j^i \frac{d^2 E_{jk}}{dT^2}}
+ {\sum_{j}^{NG} \Theta_j^i E_{jk}}
+ - \left(\frac{\partial \Lambda_k^i}{\partial T} \right)^2
+$$
+
+### Temperature-compositional cross derivative:
+
+$$
+ \frac{\partial \left(\frac{G^{E, R}}{R T} \right)}
+ {\partial n_\alpha \partial T}=
+ -\sum_k^{\mathrm{NG}} v_k^\alpha Q_k \left(\frac{\partial \Lambda_k}
+ {\partial T} - \frac{\partial \Lambda_k^\alpha}{\partial T}\right)
+ -\sum_k^{\mathrm{NG}} \left(\sum_i^{\mathrm{NC}} n_i v_k^i \right)
+ Q_k \frac{\partial^2 \Lambda_k}{\partial n_\alpha \partial T}
+$$
+
+With:
+
+$$
+ \frac{\partial^2 \Lambda_k}{\partial n_\alpha \partial T} =
+ \frac{\sum_j^{\mathrm{NG}} v_j^\alpha Q_j \frac{\partial
+ \tilde{E}_{j k}}{\partial T}}{\sum_l^{\mathrm{NC}} n_l
+ \sum_j^{\mathrm{NG}} v_j^l Q_j \tilde{E}_{j k}} -
+ \frac{\left(\sum_j^{\mathrm{NG}} v_j^\alpha Q_j \tilde{E}_{j k}\right)
+ \left(\sum_l^{\mathrm{NC}} n_l \sum_j^{\mathrm{NG}} v_j^l Q_j
+ \frac{\partial \tilde{E}_{j k}}{\partial T}\right)}
+ {\left(\sum_l^{\mathrm{NC}} n_l
+ \sum_j^{\mathrm{NG}} v_j^l Q_j \tilde{E}_{j k}\right)^2}
+$$
+
+## Examples
+### Calculating activity coefficients
+We can instantiate a [[UNIFAC]] model with a mixture ethanol-water and evaluate
+the logarithm of activity coefficients of the model for a 0.5 mole fraction of
+each, and a temperature of 298.0 K.
 
 ```fortran
 use yaeos__constants, only: pr
@@ -102,5 +308,5 @@ extension. Industrial & Engineering Chemistry Research, 30 (10), 2352–2355.
 UNIFAC Group Contribution. 6. Revision and Extension. Industrial & Engineering
 Chemistry Research, 42(1), 183–188.
 [https://doi.org/10.1021/ie020506l](https://doi.org/10.1021/ie020506l)
-9. [SINTEF - Thermopack](https://github.com/thermotools/thermopack)
+9. [SINTEF - Thermopack](https://github.com/termotools/thermopack)
 

doc/page/usage/excessmodels/uniquac.md

@@ -29,6 +29,24 @@ $$
 e_{lk}{T^2}
 $$
 
+Some of the model's terms can be simplified to reduce the complexity of the
+derivatives. Also, allows the model to be evaluated in mole vector `n` where
+some of the composition are equal to zero.
+
+$$\frac{\phi_k}{x_k} = \frac{n_T r_k}{\sum_l r_l n_l}$$
+
+$$\frac{\theta_k}{\phi_k} = \frac{q_k \sum_l r_l n_l}{r_k \sum_l q_l n_l}$$
+
+Being \(n_T \) the total number of moles in the system. The expression for the
+Excess Gibbs free energy can be rewritten as:
+
+$$ 
+\frac{G^E}{RT} = \sum_k n_k \ln \left(\frac{n_T r_k}{\sum_l r_l n_l} \right)
++ \frac{z}{2}\sum_k q_k n_k \ln \left(\frac{q_k \sum_l r_l n_l}{r_k \sum_l q_l
+n_l} \right) - \sum_k q_k n_k \ln\left(\sum_l \theta_l \tau_{lk} \right) 
+$$
+
+
 ## Temperature derivatives
 
 \(\qquad \tau_{lk}:\)
@@ -78,6 +96,8 @@ $$
 
 \(\qquad \phi_k\):
 
+\(\phi_k \) derivatives are not really needed, but we also provide them.
+
 $$
 \frac{d \phi_k}{dn_i} = \begin{cases} - \frac{{n}_{i}
 {r}_{i}^{2}}{\left(\sum_l {n}_{l} {r}_{l}\right)^{2}} +
@@ -127,48 +147,40 @@ $$
 \(\qquad G^E\):
 
 $$
-\frac{\partial \frac{G^E}{RT}}{\partial n_i} = \ln \left(\frac{\phi_i}{x_i}
-\right) + \sum_k n_k \left(\frac{\frac{d \phi_k}{dn_i}}{\phi_k} -
-\frac{\frac{dx_k}{dn_i}}{x_k}\right) +
-\frac{z}{2}{q}_{i}\ln{\left(\frac{\theta_{i}}{\phi_{i}} \right)} + \frac{z}{2}
-\sum_k {n}_{k} {q}_{k} \left(\frac{\frac{d \theta_{k}}{d {n}_{i}}}{\theta_k} -
-\frac{\frac{d \phi_{k}}{d {n}_{i}}}{\phi_k} \right) - {q}_{i}
+\frac{\partial \frac{G^E}{RT}}{\partial n_i} = \ln 
+\left(\frac{n_T r_i}{\sum_l r_l n_l} \right) + \sum_k n_k 
+\left(\frac{\sum_l r_l n_l - n_T r_i}{n_T \sum_l r_l n_l}\right) +
+\frac{z}{2} q_i \ln \left(\frac{q_i \sum_l r_l n_l}{r_i \sum_l q_l n_l} 
+\right) + \frac{z}{2} \sum_k n_k q_k \left(\frac{r_i \sum_l n_l q_l - q_i
+\sum_l n_l r_l}{(\sum_l n_l q_l) (\sum_l n_l r_l)} \right) - {q}_{i}
 \ln{\left(\sum_l \theta_{l} {\tau}_{li} \right)} - \sum_k {n}_{k} {q}_{k}
 \frac{\sum_l \frac{d \theta_{l}}{d {n}_{i}} {\tau}_{lk}}{\sum_l \theta_{l}
 {\tau}_{lk}}
 $$
 
-
 Differentiating each term of the first compositional derivative respect to
-$n_j$
+\(n_j\) we get:
 
 \(\frac{\partial \frac{G^E}{RT}}{\partial n_i \partial n_j} =\)
 
 $$
-\frac{\frac{d \phi_i}{d n_j}}{\phi_i} - \frac{\frac{d x_i}{d n_j}}{x_i} 
+\frac{\sum_l n_l r_l - n_T r_j}{n_T \sum_l n_l r_l} 
 $$
 
 $$
-+\frac{\frac{d \phi_j}{dn_i}}{\phi_j} -
-\frac{\frac{dx_j}{dn_i}}{x_j}
-+ \sum_k n_k \left(\frac{\frac{d^2\phi_k}{dn_i dn_j} \phi_k -
-  \frac{d\phi_k}{dn_i} \frac{d\phi_k}{dn_j}}{\phi_k^2} \right)
-- \sum_k n_k \left(\frac{\frac{d^2x_k}{dn_i dn_j} x_k -
-  \frac{dx_k}{dn_i} \frac{dx_k}{dn_j}}{x_k^2} \right)
++ \frac{\sum_l n_l r_l - n_T r_i}{n_T \sum_l r_l n_l}
++ \sum_k n_k \left(\frac{r_i r_j}{(\sum_l n_l r_l)^2} - \frac{1}{n_T^2} \right)
 $$
 
 $$
-+ \frac{z}{2} q_i \left( \frac{\frac{d \theta_i}{d n_j}}{\theta_i} - \frac
-{\frac{d \phi_i}{d n_j}}{\phi_i} \right)
++ \frac{z}{2} q_i \left( \frac{r_j \sum_l n_l q_l - q_j \sum_l n_l r_l}
+{(\sum_l n_l q_l)(\sum_l n_l r_l)}\right)
 $$
 
 $$
-+ \frac{z}{2} q_j \left( \frac{\frac{d \theta_j}{d n_i}}{\theta_j} - \frac
-{\frac{d \phi_j}{d n_i}}{\phi_j} \right)
-+ \frac{z}{2} \sum_k n_k q_k \left(\frac{\frac{d^2\theta_k}{dn_i dn_j} \theta_k -
-  \frac{d\theta_k}{dn_i} \frac{d\theta_k}{dn_j}}{\theta_k^2} \right)
-- \frac{z}{2} \sum_k n_k q_k \left(\frac{\frac{d^2\phi_k}{dn_i dn_j} \phi_k -
-  \frac{d\phi_k}{dn_i} \frac{d\phi_k}{dn_j}}{\phi_k^2} \right)
+\frac{z}{2} q_j \left( \frac{r_i \sum_l n_l q_l - q_i \sum_l n_l r_l}
+{(\sum_l n_l q_l)(\sum_l n_l r_l)}\right) + \frac{z}{2} \sum_k n_k q_k \left(
+\frac{q_i q_j}{(\sum_l n_l q_l)^2} - \frac{r_i r_j}{(\sum_l n_l r_l)^2} \right)
 $$
 
 $$
@@ -189,7 +201,7 @@ $$
 Example from: Gmehling et al. (2012) [2]
 
 An example of having a mixture of Water-Ethanol-Bezene at 298.15 K with 
-constant \(\Delta U\) [K]:
+constant \(\frac{\Delta U}{R}\) [K]:
 
 |Water|Ethanol|Benzene|
 |---------|--------|---------|