Difference between revisions of "DOW Experimental Design"
RobertLampel (Talk | contribs) (→Chemical background) |
RobertLampel (Talk | contribs) (→Chemical background) |
||
(19 intermediate revisions by 2 users not shown) | |||
Line 47: | Line 47: | ||
</math> | </math> | ||
</p> | </p> | ||
− | Here <math>[ \ ]</math> denotes the concentration of the species in <math> | + | Here <math>[ \ ]</math> denotes the concentration of the species in <math>mol/kg</math>. |
By assuming the rapid acid-base reactions are at equilibrium, the equilibrium constants <math>K_1,K_2,K_3</math> can be defined as | By assuming the rapid acid-base reactions are at equilibrium, the equilibrium constants <math>K_1,K_2,K_3</math> can be defined as | ||
<p> | <p> | ||
Line 85: | Line 85: | ||
<math> | <math> | ||
\begin{align} | \begin{align} | ||
− | \frac{d[MBMH]}{dt} &= k_1 \left[ M^- \right] \left[ BM \right] | + | \frac{d[MBMH]}{dt} &= k_1 \left[ M^- \right] \left[ BM \right] - k_{-1} \left[ MBM^- \right] && \quad (g)\\ |
\frac{d[HA]}{dt} &= k_2 \left[ A^- \right] \left[ BM \right] && \quad (h)\\ | \frac{d[HA]}{dt} &= k_2 \left[ A^- \right] \left[ BM \right] && \quad (h)\\ | ||
\frac{d[HABM]}{dt} &= k_2 \left[ A^- \right] \left[ BM \right] + k_3 \left[ M^- \right] \left[ AB \right] - k_{-3} \left[ ABM^- \right] && \quad (i) | \frac{d[HABM]}{dt} &= k_2 \left[ A^- \right] \left[ BM \right] + k_3 \left[ M^- \right] \left[ AB \right] - k_{-3} \left[ ABM^- \right] && \quad (i) | ||
Line 110: | Line 110: | ||
</math> | </math> | ||
</p> | </p> | ||
− | Here <math>R \approx | + | Here <math>R \approx 1.98720425864083 \ cal/(K \cdot mol) </math> is the gas constant and <math>T</math> is the reaction temperature in Kelvins. The parameters <math>\alpha_i</math>, given in <math>mol/( kg \cdot h)</math>, represent the pre-exponential factors and the <math>E_i</math>, with unit <math>cal/mol</math>, are the activation energies. |
== Mathematical formulation == | == Mathematical formulation == | ||
Line 120: | Line 120: | ||
\dot{y}_1 &= -k_2 y_8 y_2 && \quad (1),(h) \\ | \dot{y}_1 &= -k_2 y_8 y_2 && \quad (1),(h) \\ | ||
\dot{y}_2 &= -k_1 y_6 y_2 + k_{-1} y_{10} - k_2 y_8 y_2 && \quad (2),(e) \\ | \dot{y}_2 &= -k_1 y_6 y_2 + k_{-1} y_{10} - k_2 y_8 y_2 && \quad (2),(e) \\ | ||
− | \dot{y}_3 &= | + | \dot{y}_3 &= k_2 y_8 y_2 + k_1 y_6 y_4 - \frac{1}{2} k_{-1} y_9 && \quad (3),(i) \\ |
\dot{y}_4 &= -k_1 y_6 y_4 + \frac{1}{2} k_{-1} y_9 && \quad (4),(f) \\ | \dot{y}_4 &= -k_1 y_6 y_4 + \frac{1}{2} k_{-1} y_9 && \quad (4),(f) \\ | ||
− | \dot{y}_5 &= | + | \dot{y}_5 &= k_1 y_6 y_2 - k_{-1} y_{10} && \quad (5),(g) \\ |
\dot{y}_6 &= -k_1 (y_6 y_2 + y_6 y_4) + k_{-1} (y_{10} + \frac{1}{2} y_9) && \quad (6),(d) \\ | \dot{y}_6 &= -k_1 (y_6 y_2 + y_6 y_4) + k_{-1} (y_{10} + \frac{1}{2} y_9) && \quad (6),(d) \\ | ||
y_7 &= -\left[ Q^+ \right] + y_6 + y_8 + y_9 + y_{10} && \quad (7),(j)\\ | y_7 &= -\left[ Q^+ \right] + y_6 + y_8 + y_9 + y_{10} && \quad (7),(j)\\ | ||
Line 131: | Line 131: | ||
</math> | </math> | ||
</p> | </p> | ||
− | Here the letters in parentheses stand for the corresponding chemical process and the quantity <math>\left[Q^+\right]=0. | + | Here the letters in parentheses stand for the corresponding chemical process and the quantity <math>\left[Q^+\right]=0.0131</math> is a constant during the reaction. |
The nine parameters form the vector | The nine parameters form the vector | ||
<p> | <p> | ||
Line 174: | Line 174: | ||
<tr> | <tr> | ||
<td style="text-align: center; padding:5pt"> <math> \alpha_1 </math> </td> | <td style="text-align: center; padding:5pt"> <math> \alpha_1 </math> </td> | ||
− | <td style="text-align: center; padding:5pt"> <math>2.0 \ | + | <td style="text-align: center; padding:5pt"> <math>2.0 \cdot 10^{13} \ mol / (kg \cdot h)</math> </td> |
</tr> | </tr> | ||
<tr> | <tr> | ||
− | <td style="text-align: center; padding:5pt"><math> | + | <td style="text-align: center; padding:5pt"><math>\alpha_2</math></td> |
− | <td style="text-align: center; padding:5pt"><math>2.0 \ | + | <td style="text-align: center; padding:5pt"><math>2.0 \cdot 10^{13} \ mol / (kg \cdot h)</math></td> |
</tr> | </tr> | ||
<tr> | <tr> | ||
− | <td style="text-align: center; padding:5pt"><math>\ | + | <td style="text-align: center; padding:5pt"><math>\alpha_{-1}</math></td> |
− | <td style="text-align: center; padding:5pt"><math> | + | <td style="text-align: center; padding:5pt"><math>4.3 \cdot 10^{15} \ mol / (kg \cdot h)</math></td> |
</tr> | </tr> | ||
<tr> | <tr> | ||
− | <td style="text-align: center; padding:5pt"><math> | + | <td style="text-align: center; padding:5pt"><math>E_1</math></td> |
− | <td style="text-align: center; padding:5pt"><math>2.0 \ | + | <td style="text-align: center; padding:5pt"><math>2.0 \cdot 10^{4} \ cal / mol</math></td> |
</tr> | </tr> | ||
<tr> | <tr> | ||
− | <td style="text-align: center; padding:5pt"><math> | + | <td style="text-align: center; padding:5pt"><math>E_2</math></td> |
− | <td style="text-align: center; padding:5pt"><math> | + | <td style="text-align: center; padding:5pt"><math>2.0 \cdot 10^{4} \ cal / mol </math></td> |
</tr> | </tr> | ||
<tr> | <tr> | ||
<td style="text-align: center; padding:5pt"><math>E_{-1}</math></td> | <td style="text-align: center; padding:5pt"><math>E_{-1}</math></td> | ||
− | <td style="text-align: center; padding:5pt"><math>2.0 \ | + | <td style="text-align: center; padding:5pt"><math>2.0 \cdot 10^{4} \ cal / mol</math></td> |
</tr> | </tr> | ||
<tr> | <tr> | ||
<td style="text-align: center; padding:5pt"><math>K_1</math></td> | <td style="text-align: center; padding:5pt"><math>K_1</math></td> | ||
− | <td style="text-align: center; padding:5pt"><math>1.0\cdot 10^{-17}</math></td> | + | <td style="text-align: center; padding:5pt"><math>1.0\cdot 10^{-17} \ mol / kg</math></td> |
</tr> | </tr> | ||
<tr> | <tr> | ||
<td style="text-align: center; padding:5pt"><math>K_2</math></td> | <td style="text-align: center; padding:5pt"><math>K_2</math></td> | ||
− | <td style="text-align: center; padding:5pt"><math>1.0\cdot 10^{-11}</math></td> | + | <td style="text-align: center; padding:5pt"><math>1.0\cdot 10^{-11} \ mol / kg</math></td> |
</tr> | </tr> | ||
<tr> | <tr> | ||
<td style="text-align: center; padding:5pt"><math>K_3</math></td> | <td style="text-align: center; padding:5pt"><math>K_3</math></td> | ||
− | <td style="text-align: center; padding:5pt"><math>1.0\cdot 10^{-17}</math></td> | + | <td style="text-align: center; padding:5pt"><math>1.0\cdot 10^{-17} \ mol / kg</math></td> |
</tr> | </tr> | ||
</table> | </table> | ||
− | + | Note that for the calculations all temperatures given in <math>^{\circ}C</math> have to be rescaled to <math>K</math> by adding <math>273.15</math>. | |
There are three datasets for different temperatures <math>T</math>, with corresponding starting values | There are three datasets for different temperatures <math>T</math>, with corresponding starting values | ||
Line 263: | Line 263: | ||
\begin{array}{l} | \begin{array}{l} | ||
y_5 = 0 \\ | y_5 = 0 \\ | ||
− | y_6 = 0. | + | y_6 = 0.00131 \\ |
y_7 = \frac{1}{2} \cdot \left( -K_2 + \sqrt{K_2^2 + 4K_2 y_1(0)} \right) \\ | y_7 = \frac{1}{2} \cdot \left( -K_2 + \sqrt{K_2^2 + 4K_2 y_1(0)} \right) \\ | ||
y_8 = y_7 \\ | y_8 = y_7 \\ | ||
Line 271: | Line 271: | ||
</math> | </math> | ||
</p> | </p> | ||
+ | |||
+ | To reduce the intercorrelation between the parameters in the rate constants, we apply the following reparametrization (cf. [[#Kristensen2004|[4]]].): | ||
+ | <p> | ||
+ | <math> | ||
+ | \begin{align} | ||
+ | k_i &= \alpha_i \cdot \exp \left( - \frac{E_i}{RT} \right) \\ | ||
+ | &= k_{0,i} \cdot \exp \left( - \frac{E_i}{R} \cdot \left( \frac{1}{T} - \frac{1}{T_0} \right) \right), \quad i=1,2,-1 | ||
+ | \end{align} | ||
+ | </math> | ||
+ | </p> | ||
+ | in which <math>k_{0,i} = \alpha_i \cdot \exp\left( - \frac{E_i}{RT_0} \right)</math>. The reference temperature in <math>T_0</math> is chosen as the average over all performed experiments, i.e., <math>T_0 = 69^\circ C</math>. Additionally, we add a logarithmic transformation, which gives rise to the following transformed starting values: | ||
+ | |||
+ | <table style="border-collapse: collapse;" border="1."> | ||
+ | <tr> | ||
+ | <td style="text-align: center; padding:5pt"> <math> \ln k_{0,1} </math> </td> | ||
+ | <td style="text-align: center; padding:5pt"> <math>1.194</math> </td> | ||
+ | </tr> | ||
+ | <tr> | ||
+ | <td style="text-align: center; padding:5pt"><math>\ln k_{0,2}</math></td> | ||
+ | <td style="text-align: center; padding:5pt"><math>1.194</math></td> | ||
+ | </tr> | ||
+ | <tr> | ||
+ | <td style="text-align: center; padding:5pt"><math>\ln k_{0,-1}</math></td> | ||
+ | <td style="text-align: center; padding:5pt"><math>6.565</math></td> | ||
+ | </tr> | ||
+ | <tr> | ||
+ | <td style="text-align: center; padding:5pt"><math>E_1</math></td> | ||
+ | <td style="text-align: center; padding:5pt"><math>2.0 \cdot 10^{4}</math></td> | ||
+ | </tr> | ||
+ | <tr> | ||
+ | <td style="text-align: center; padding:5pt"><math>E_2</math></td> | ||
+ | <td style="text-align: center; padding:5pt"><math>2.0 \cdot 10^{4} </math></td> | ||
+ | </tr> | ||
+ | <tr> | ||
+ | <td style="text-align: center; padding:5pt"><math>E_{-1}</math></td> | ||
+ | <td style="text-align: center; padding:5pt"><math>2.0 \cdot 10^{4}</math></td> | ||
+ | </tr> | ||
+ | <tr> | ||
+ | <td style="text-align: center; padding:5pt"><math>\ln K_1</math></td> | ||
+ | <td style="text-align: center; padding:5pt"><math>-34.54</math></td> | ||
+ | </tr> | ||
+ | <tr> | ||
+ | <td style="text-align: center; padding:5pt"><math>\ln K_2</math></td> | ||
+ | <td style="text-align: center; padding:5pt"><math>-25.33</math></td> | ||
+ | </tr> | ||
+ | <tr> | ||
+ | <td style="text-align: center; padding:5pt"><math>\ln K_3</math></td> | ||
+ | <td style="text-align: center; padding:5pt"><math>-39.14</math></td> | ||
+ | </tr> | ||
+ | </table> | ||
== Optimal Experimental Design Problem == | == Optimal Experimental Design Problem == | ||
Line 365: | Line 415: | ||
<span id="OEDUDE">[2]</span> "Optimal Experimental Design for Universal Differential Equations" by C. Plate, C.J. Martensen and S. Sager <br> | <span id="OEDUDE">[2]</span> "Optimal Experimental Design for Universal Differential Equations" by C. Plate, C.J. Martensen and S. Sager <br> | ||
<span id="Stortelder1998">[3]</span> "Parameter estimation in nonlinear systems" by W.J.H. Stortelder <br> | <span id="Stortelder1998">[3]</span> "Parameter estimation in nonlinear systems" by W.J.H. Stortelder <br> | ||
+ | <span id="Kristensen2004">[4]</span> "Parameter Estimation in Nonlinear Dynamical Systems" by Morten Rode Kristensen <br> | ||
<!--List of all categories this page is part of. List characterization of solution behavior, model properties, ore presence of implementation details (e.g., AMPL for AMPL model) here --> | <!--List of all categories this page is part of. List characterization of solution behavior, model properties, ore presence of implementation details (e.g., AMPL for AMPL model) here --> |
Latest revision as of 08:34, 18 December 2024
DOW Experimental Design | |
---|---|
State dimension: | 1 |
Differential states: | 11 |
Discrete control functions: | 2 |
Path constraints: | 4 |
Interior point equalities: | 11 |
The DOW Experimental Design problem models the OED problem for the parameter estimation problem formulated by the DOW Chemical Co. in 1981. The following formulation is taken from [1].
Contents
Chemical background
The chemical species are disguised for proprietary reasons and the desired reaction is given by , where is the desired product. The reactions are described as follows:
Slow Kinetic Reactions:
Acid-Base Reactions:
In order to devise a model to account for these reactions, it is first necessary to distinguish between the overall concentration of a species and the concentration of its neutral form. Overall concentrations are defined for three components based on neutral and ionic species
Here denotes the concentration of the species in . By assuming the rapid acid-base reactions are at equilibrium, the equilibrium constants can be defined as
The anionic species may then be represented by
Material balance equations for the three reactants in the slow kinetic reactions yield:
From stoichiometry, rate expressions can also be written for the total species:
An electroneutrality constraint gives the hydrogen ion concentration as
Based on similarities of reacting species, we assume
With these assumptions, the three rate constants and must be estimated. Each of these can be fitted with two adjustable model parameters, assuming an Arrhenius temperature dependence. That is
Here is the gas constant and is the reaction temperature in Kelvins. The parameters , given in , represent the pre-exponential factors and the , with unit , are the activation energies.
Mathematical formulation
The chemical processes can be expressed mathematically as six differential equations and four algebraic equations:
Here the letters in parentheses stand for the corresponding chemical process and the quantity is a constant during the reaction. The nine parameters form the vector
The predicted concentrations form the vector
Let denote the right hand side of equation for . We reformulate the last four algebraic equations as differential ones:
The right hand sides of and are summarized as the vector-valued function . Moreover, let
Parameters
The initial parameter estimates are:
Note that for the calculations all temperatures given in have to be rescaled to by adding .
There are three datasets for different temperatures , with corresponding starting values
The initial model conditions in addition to those given in the data sets are:
To reduce the intercorrelation between the parameters in the rate constants, we apply the following reparametrization (cf. [4].):
in which . The reference temperature in is chosen as the average over all performed experiments, i.e., . Additionally, we add a logarithmic transformation, which gives rise to the following transformed starting values:
Optimal Experimental Design Problem
To be specified.
We are interested in when to measure (with an upper bound on the measuring time for each observable). We define
In this approach, we add the so-called sensitivities . For the differential equations this means
Now we formulate the OED problem as described in [2].
Here is the observed function. The evolution of the symmetric matrix is given by the weighted sum of observability Gramians for each observed function of states. The weights are the (binary) sampling decisions, where denotes the decision to perform a measurement at time .
Miscellaneous and Further Reading
To be specified.
References
[1] "Nonlinear Parameter Estimation: a Case Study Comparison" by L. T. Biegler and J. J. Damiano
[2] "Optimal Experimental Design for Universal Differential Equations" by C. Plate, C.J. Martensen and S. Sager
[3] "Parameter estimation in nonlinear systems" by W.J.H. Stortelder
[4] "Parameter Estimation in Nonlinear Dynamical Systems" by Morten Rode Kristensen