Show code cell content
###############################################################################
# The Institute for the Design of Advanced Energy Systems Integrated Platform
# Framework (IDAES IP) was produced under the DOE Institute for the
# Design of Advanced Energy Systems (IDAES).
#
# Copyright (c) 2018-2026 by the software owners: The Regents of the
# University of California, through Lawrence Berkeley National Laboratory,
# National Technology & Engineering Solutions of Sandia, LLC, Carnegie Mellon
# University, West Virginia University Research Corporation, et al.
# All rights reserved. Please see the files COPYRIGHT.md and LICENSE.md
# for full copyright and license information.
###############################################################################
Parameter Estimation Using Flash Unit Model#
Author: Jaffer Ghouse
Maintainer: Stephen Cini
Updated: 2026-06-11
In this module, we will be using Pyomo’s parmest tool in conjunction with IDAES models for parameter estimation. We demonstrate these tools by estimating the parameters associated with the NRTL property model for a benzene-toluene mixture. The NRTL model has 2 sets of parameters: the non-randomness parameter (alpha_ij) and the binary interaction parameter (tau_ij), where i and j is the pure component species. In this example, we will be only estimate the binary interaction parameter (tau_ij) for a given dataset. When estimating parameters associated with the property package, IDAES provides the flexibility of doing the parameter estimation by just using the state block or by using a unit model with a specified property package. This module will demonstrate parameter estimation by using the flash unit model with the NRTL property package.
We will complete the following tasks:
Set up a method to return an initialized model
Set up the parameter estimation problem using
parmestAnalyze the results
Demonstrate advanced features from
parmest
Key links to documentation:#
# Todo: import ConcreteModel, value, and Suffix from pyomo.environ
from pyomo.environ import ConcreteModel, value, Suffix
# Todo: import FlowsheetBlock from idaes.core
from idaes.core import FlowsheetBlock
# Todo: import Flash unit model from idaes.models.unit_models
from idaes.models.unit_models import Flash
In the next cell, we will be importing the parameter block that we will be using in this module and the idaes logger.
from idaes.models.properties.activity_coeff_models.BTX_activity_coeff_VLE import (
BTXParameterBlock,
)
import idaes.logger as idaeslog
In the next cell, we import parmest from Pyomo and the pandas package. We need pandas as parmest uses pandas.dataframe for handling the input data and the results.
import pyomo.contrib.parmest.parmest as parmest
import pandas as pd
Setting up an initialized model#
We need to provide a method that returns an initialized model to the parmest tool in Pyomo.
def NRTL_model(data):
# Todo: Create a ConcreteModel object
m = ConcreteModel()
# Todo: Create FlowsheetBlock object
m.fs = FlowsheetBlock(dynamic=False)
# Todo: Create a properties parameter object with the following options:
# "valid_phase": ('Liq', 'Vap')
# "activity_coeff_model": 'NRTL'
m.fs.properties = BTXParameterBlock(
valid_phase=("Liq", "Vap"), activity_coeff_model="NRTL"
)
m.fs.flash = Flash(property_package=m.fs.properties)
# Initialize at a certain inlet condition
m.fs.flash.inlet.flow_mol.fix(1)
m.fs.flash.inlet.temperature.fix(368)
m.fs.flash.inlet.pressure.fix(101325)
m.fs.flash.inlet.mole_frac_comp[0, "benzene"].fix(0.5)
m.fs.flash.inlet.mole_frac_comp[0, "toluene"].fix(0.5)
# Set Flash unit specifications
m.fs.flash.heat_duty.fix(0)
m.fs.flash.deltaP.fix(0)
# Fix NRTL specific variables
# alpha values (set at 0.3)
m.fs.properties.alpha["benzene", "benzene"].fix(0)
m.fs.properties.alpha["benzene", "toluene"].fix(0.3)
m.fs.properties.alpha["toluene", "toluene"].fix(0)
m.fs.properties.alpha["toluene", "benzene"].fix(0.3)
# initial tau values
m.fs.properties.tau["benzene", "benzene"].fix(0)
m.fs.properties.tau["benzene", "toluene"].fix(-0.9)
m.fs.properties.tau["toluene", "toluene"].fix(0)
m.fs.properties.tau["toluene", "benzene"].fix(1.4)
# Initialize the flash unit
m.fs.flash.initialize(outlvl=idaeslog.INFO_LOW)
# Fix at actual temperature
if isinstance(data, dict) or isinstance(data, pd.Series):
m.fs.flash.inlet.temperature.fix(float(data["temperature"]))
elif isinstance(data, pd.DataFrame):
m.fs.flash.inlet.temperature.fix(float(data.iloc[0]["temperature"]))
else:
raise ValueError("Unrecognized data type.")
# Set bounds on variables to be estimated
m.fs.properties.tau["benzene", "toluene"].setlb(-5)
m.fs.properties.tau["benzene", "toluene"].setub(5)
m.fs.properties.tau["toluene", "benzene"].setlb(-5)
m.fs.properties.tau["toluene", "benzene"].setub(5)
# Return initialized flash model
return m
Parameter estimation using parmest#
In addition to providing a method to return an initialized model, the parmest tool needs the following:
Experiment class to set up and label model with suffixes
Dataset with multiple scenarios - organized into an experiment list
Here we build an experiment class to label our model problem for parameter estimation. The labels are defined as a Suffix, and the main labels for our model are experiment_outputs, unknown_parameters, and measurement_error.
For this problem, the error will be computed for the mole fraction of benzene in the vapor and liquid phase between the model prediction and data. The experimental_outputs will therefore be the mole fraction of benzene in the two phases.
In this example, we only estimate the binary interaction parameter (tau_ij). Given that this variable is usually indexed as tau_ij = Var(component_list, component_list), there are 2*2=4 degrees of freedom. However, when i=j, the binary interaction parameter is 0. Therefore, in this problem, we estimate the binary interaction parameter for the following variables only:
fs.properties.tau[‘benzene’, ‘toluene’]
fs.properties.tau[‘toluene’, ‘benzene’]
As shown below, these model components are used as our unknown_parameters.
We define measurement_error as none so parmest calculates the value internally based on the experimental outputs. Refer to https://pyomo.readthedocs.io/en/stable/explanation/analysis/parmest/driver.html for more information.
# Build an experiment class to take advantage of new parmest interface
from pyomo.contrib.parmest.experiment import Experiment
class NRTLExperiment(Experiment):
"""Experiment class for parameter estimation of NRTL model using parmest"""
def __init__(self, data, meas_error=None):
"""Initialize the NRTLExperiment class
Args:
data: DataFrame containing the experimental data
meas_error: Measurement error for the data (optional)
"""
self.model = None
self.data = data
self.meas_error = meas_error
def create_model(self):
"""Create the Pyomo model for the NRTL parameter estimation problem"""
self.model = NRTL_model(self.data)
def label_model(self):
m = self.model
# Add experiment outputs to the model for easier access
m.experiment_outputs = Suffix(direction=Suffix.LOCAL)
m.experiment_outputs.update(
[
(
m.fs.flash.liq_outlet.mole_frac_comp[0, "benzene"],
self.data["liq_benzene"],
),
(
m.fs.flash.vap_outlet.mole_frac_comp[0, "benzene"],
self.data["vap_benzene"],
),
]
)
m.measurement_error = Suffix(direction=Suffix.LOCAL)
m.measurement_error.update(
[
(m.fs.flash.liq_outlet.mole_frac_comp[0, "benzene"], self.meas_error),
(m.fs.flash.vap_outlet.mole_frac_comp[0, "benzene"], self.meas_error),
]
)
# Add unknown parameters to the model for easier access
m.unknown_parameters = Suffix(direction=Suffix.LOCAL)
m.unknown_parameters.update(
(k, value(k))
for k in [
m.fs.properties.tau["benzene", "toluene"],
m.fs.properties.tau["toluene", "benzene"],
]
)
def get_labeled_model(self):
"""Return the labeled model"""
if self.model is None:
self.create_model()
self.label_model()
return self.model
Pyomo’s parmest tool supports the following data formats:
pandas dataframe
list of dictionaries
list of json file names.
Please see the documentation for more details.
For this example, we load data from the csv file BT_NRTL_dataset.csv. The dataset consists of fifty data points which provide the mole fraction of benzene in the vapor and liquid phase as a function of temperature.
# Load data from csv
data = pd.read_csv("BT_NRTL_dataset.csv")
# Display the dataset
display(data)
| temperature | liq_benzene | vap_benzene | |
|---|---|---|---|
| 0 | 365.500000 | 0.480953 | 0.692110 |
| 1 | 365.617647 | 0.462444 | 0.667699 |
| 2 | 365.735294 | 0.477984 | 0.692441 |
| 3 | 365.852941 | 0.440547 | 0.640336 |
| 4 | 365.970588 | 0.427421 | 0.623328 |
| 5 | 366.088235 | 0.442725 | 0.647796 |
| 6 | 366.205882 | 0.434374 | 0.637691 |
| 7 | 366.323529 | 0.444642 | 0.654933 |
| 8 | 366.441176 | 0.427132 | 0.631229 |
| 9 | 366.558824 | 0.446301 | 0.661743 |
| 10 | 366.676471 | 0.438004 | 0.651591 |
| 11 | 366.794118 | 0.425320 | 0.634814 |
| 12 | 366.911765 | 0.439435 | 0.658047 |
| 13 | 367.029412 | 0.435655 | 0.654539 |
| 14 | 367.147059 | 0.401350 | 0.604987 |
| 15 | 367.264706 | 0.397862 | 0.601703 |
| 16 | 367.382353 | 0.415821 | 0.630930 |
| 17 | 367.500000 | 0.420667 | 0.640380 |
| 18 | 367.617647 | 0.391683 | 0.598214 |
| 19 | 367.735294 | 0.404903 | 0.620432 |
| 20 | 367.852941 | 0.409563 | 0.629626 |
| 21 | 367.970588 | 0.389488 | 0.600722 |
| 22 | 368.000000 | 0.396789 | 0.612483 |
| 23 | 368.088235 | 0.398162 | 0.616106 |
| 24 | 368.205882 | 0.362340 | 0.562505 |
| 25 | 368.323529 | 0.386958 | 0.602680 |
| 26 | 368.441176 | 0.363643 | 0.568210 |
| 27 | 368.558824 | 0.368118 | 0.577072 |
| 28 | 368.676471 | 0.384098 | 0.604078 |
| 29 | 368.794118 | 0.353605 | 0.557925 |
| 30 | 368.911765 | 0.346474 | 0.548445 |
| 31 | 369.029412 | 0.350741 | 0.556996 |
| 32 | 369.147059 | 0.362347 | 0.577286 |
| 33 | 369.264706 | 0.362578 | 0.579519 |
| 34 | 369.382353 | 0.340765 | 0.546411 |
| 35 | 369.500000 | 0.337462 | 0.542857 |
| 36 | 369.617647 | 0.355729 | 0.574083 |
| 37 | 369.735294 | 0.348679 | 0.564513 |
| 38 | 369.852941 | 0.338187 | 0.549284 |
| 39 | 369.970588 | 0.324360 | 0.528514 |
| 40 | 370.088235 | 0.310753 | 0.507964 |
| 41 | 370.205882 | 0.311037 | 0.510055 |
| 42 | 370.323529 | 0.311263 | 0.512055 |
| 43 | 370.441176 | 0.308081 | 0.508437 |
| 44 | 370.558824 | 0.308224 | 0.510293 |
| 45 | 370.676471 | 0.318148 | 0.528399 |
| 46 | 370.794118 | 0.308334 | 0.513728 |
| 47 | 370.911765 | 0.317937 | 0.531410 |
| 48 | 371.029412 | 0.289149 | 0.484824 |
| 49 | 371.147059 | 0.298637 | 0.502318 |
We define the exp_list by splitting the data into individual experiments, or data points.
# Loop through the dataset and create an experiment for each row of data
exp_list = []
for i in range(data.shape[0]):
exp_list.append(NRTLExperiment(data.iloc[i]))
We are now ready to set up the parameter estimation problem. We will create a parameter estimation object called pest. As shown below, we pass the experiment list, and an objective function to the Estimator method. tee=True will print the solver output after solving the parameter estimation problem.
import logging
idaeslog.getIdaesLogger("core.property_meta").setLevel(logging.ERROR)
pest = parmest.Estimator(exp_list, obj_function="SSE", tee=True)
obj_value, parameters = pest.theta_est()
Ipopt 3.13.2:
******************************************************************************
This program contains Ipopt, a library for large-scale nonlinear optimization.
Ipopt is released as open source code under the Eclipse Public License (EPL).
For more information visit http://projects.coin-or.org/Ipopt
This version of Ipopt was compiled from source code available at
https://github.com/IDAES/Ipopt as part of the Institute for the Design of
Advanced Energy Systems Process Systems Engineering Framework (IDAES PSE
Framework) Copyright (c) 2018-2019. See https://github.com/IDAES/idaes-pse.
This version of Ipopt was compiled using HSL, a collection of Fortran codes
for large-scale scientific computation. All technical papers, sales and
publicity material resulting from use of the HSL codes within IPOPT must
contain the following acknowledgement:
HSL, a collection of Fortran codes for large-scale scientific
computation. See http://www.hsl.rl.ac.uk.
******************************************************************************
This is Ipopt version 3.13.2, running with linear solver ma27.
Number of nonzeros in equality constraint Jacobian...: 10950
Number of nonzeros in inequality constraint Jacobian.: 0
Number of nonzeros in Lagrangian Hessian.............: 6600
Total number of variables............................: 2952
variables with only lower bounds: 150
variables with lower and upper bounds: 600
variables with only upper bounds: 0
Total number of equality constraints.................: 2950
Total number of inequality constraints...............: 0
inequality constraints with only lower bounds: 0
inequality constraints with lower and upper bounds: 0
inequality constraints with only upper bounds: 0
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
0 6.0671019e-03 5.63e+02 1.20e-08 -1.0 0.00e+00 - 0.00e+00 0.00e+00 0
1 9.0714942e-04 1.37e+03 1.61e+01 -1.0 1.37e+04 - 9.82e-01 1.00e+00h 1
2 1.3219937e-03 1.17e+04 1.62e+01 -1.0 5.23e+03 - 3.41e-01 1.65e-01h 3
3 1.3161671e-03 1.11e+04 4.02e+01 -1.0 3.96e+02 -4.0 8.19e-01 1.25e-01h 4
4 1.3154170e-03 1.11e+04 4.32e+01 -1.0 3.47e+02 -4.5 9.90e-01 4.43e-02h 5
5 9.4424343e-04 1.04e+04 4.05e+01 -1.0 1.16e+04 - 9.33e-01 5.61e-02h 5
6 9.4571258e-04 1.11e+04 4.85e+01 -1.0 3.13e+02 -5.0 8.74e-01 1.25e-01h 4
7 9.4862862e-04 1.11e+04 4.84e+01 -1.0 2.41e+03 -5.4 1.72e-01 1.50e-03h 8
8 9.5448936e-04 1.10e+04 4.81e+01 -1.0 4.11e+02 -5.0 1.00e+00 1.97e-02h 6
9 9.5650047e-04 1.10e+04 4.81e+01 -1.0 5.14e+02 -4.6 3.46e-01 5.24e-03h 7
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
10 9.5730419e-04 1.10e+04 4.81e+01 -1.0 9.78e+02 -4.2 4.40e-01 1.08e-03h 8
11 9.5774182e-04 1.10e+04 4.84e+01 -1.0 2.76e+02 -3.7 2.91e-01 2.99e-03h 7
12 9.4292855e-04 7.77e+04 4.09e+04 -1.0 5.38e+02 -2.4 5.48e-02 6.44e-02w 1
13 2.4318219e-01 2.38e+07 4.22e+14 -1.0 1.28e+06 - 2.55e-02 3.49e-02w 1
14 2.8902019e-01 1.64e+07 2.92e+14 -1.0 1.11e+05 -2.9 9.54e-01 3.01e-01w 1
15 9.5768293e-04 1.10e+04 4.85e+01 -1.0 2.80e+05 -3.4 5.48e-02 2.52e-04h 8
16 9.5769345e-04 1.10e+04 4.85e+01 -1.0 4.33e+02 -2.9 7.95e-02 1.91e-04h 9
17 9.5763809e-04 1.10e+04 4.85e+01 -1.0 6.35e+02 -2.5 4.94e-02 2.08e-04h 9
18 9.5734836e-04 1.10e+04 4.87e+01 -1.0 3.22e+02 -3.0 1.00e+00 1.22e-03h 8
19 9.3628114e-04 1.09e+04 4.62e+01 -1.0 4.22e+02 -3.5 5.16e-01 5.87e-02h 5
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
20 7.1894543e-04 4.02e+03 4.95e+02 -1.0 3.49e+02 -3.9 1.00e+00 1.00e+00h 1
21 6.7062834e-04 1.58e+02 2.21e+02 -1.0 1.78e+02 -4.4 1.00e+00 1.00e+00h 1
22 6.6635123e-04 4.17e+01 1.12e+01 -1.0 2.44e+02 -4.9 1.00e+00 1.00e+00h 1
23 6.6729070e-04 1.43e+00 9.69e-01 -1.0 1.08e+01 -5.4 1.00e+00 1.00e+00h 1
24 6.6785681e-04 3.58e-01 5.25e-03 -1.7 6.65e+00 -5.8 1.00e+00 1.00e+00h 1
25 6.2552355e-04 9.80e+02 6.48e-02 -3.8 7.92e+02 - 9.22e-01 1.00e+00h 1
26 5.9105707e-04 8.13e+00 6.13e-04 -3.8 5.66e+02 - 1.00e+00 1.00e+00h 1
27 6.1313801e-04 6.24e+01 8.08e-06 -3.8 1.70e+02 - 1.00e+00 1.00e+00h 1
28 6.1208051e-04 3.58e-01 1.18e-08 -3.8 1.88e+00 - 1.00e+00 1.00e+00h 1
29 5.9010560e-04 1.58e+01 1.00e-05 -5.7 1.08e+02 - 1.00e+00 1.00e+00h 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
30 5.9008319e-04 2.74e-02 7.90e-07 -5.7 7.54e-02 -6.3 1.00e+00 1.00e+00h 1
31 5.3570104e-04 1.48e+04 8.70e-04 -5.7 1.81e+02 - 1.00e+00 1.00e+00h 1
32 5.1618763e-04 6.29e+02 8.06e-05 -5.7 2.55e+02 - 1.00e+00 1.00e+00h 1
33 5.5250565e-04 2.88e+01 2.58e-04 -5.7 1.81e+04 - 6.26e-01 6.25e-02h 5
34 5.2209909e-04 2.88e+02 2.02e-04 -5.7 1.89e+03 - 1.00e+00 1.00e+00h 1
35 5.0707798e-04 5.60e+01 2.02e-04 -5.7 2.93e+03 - 1.00e+00 1.00e+00h 1
36 5.0765648e-04 7.87e+00 1.91e-05 -5.7 1.02e+02 - 1.00e+00 1.00e+00h 1
37 5.0740606e-04 2.57e+00 1.88e-05 -5.7 6.96e+01 - 1.00e+00 1.00e+00h 1
38 5.0764164e-04 2.55e+00 1.76e-05 -5.7 6.81e+01 - 1.00e+00 1.00e+00h 1
39 5.0740672e-04 2.55e+00 1.88e-05 -5.7 6.80e+01 - 1.00e+00 1.00e+00h 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
40 5.0764162e-04 2.54e+00 1.76e-05 -5.7 6.80e+01 - 1.00e+00 1.00e+00h 1
41 5.0740672e-04 2.55e+00 1.88e-05 -5.7 6.80e+01 - 1.00e+00 1.00e+00h 1
42 5.0782847e-04 1.43e-01 1.83e-05 -5.7 6.80e+01 - 1.00e+00 1.00e+00H 1
43 5.0751252e-04 1.58e-02 1.93e-05 -5.7 6.88e+01 - 1.00e+00 1.00e+00H 1
44 5.0783218e-04 4.79e-03 1.85e-05 -5.7 6.92e+01 - 1.00e+00 1.00e+00H 1
45 5.0751252e-04 1.60e-02 1.93e-05 -5.7 6.89e+01 - 1.00e+00 1.00e+00H 1
46 5.0783222e-04 1.23e-04 1.85e-05 -5.7 6.92e+01 - 1.00e+00 1.00e+00H 1
47 5.0751252e-04 1.60e-02 1.93e-05 -5.7 6.89e+01 - 1.00e+00 1.00e+00H 1
48 5.0783216e-04 1.24e-04 1.85e-05 -5.7 6.92e+01 - 1.00e+00 1.00e+00H 1
49 5.0751765e-04 3.20e-03 1.93e-05 -5.7 6.89e+01 - 1.00e+00 1.00e+00H 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
50 5.0764151e-04 2.56e+00 1.77e-05 -5.7 6.93e+01 - 1.00e+00 1.00e+00h 1
51 5.0740671e-04 2.55e+00 1.88e-05 -5.7 6.81e+01 - 1.00e+00 1.00e+00h 1
52 5.0764162e-04 2.54e+00 1.76e-05 -5.7 6.80e+01 - 1.00e+00 1.00e+00h 1
53 5.0740672e-04 2.55e+00 1.88e-05 -5.7 6.80e+01 - 1.00e+00 1.00e+00h 1
54 5.0764162e-04 2.54e+00 1.76e-05 -5.7 6.80e+01 - 1.00e+00 1.00e+00h 1
55 5.0748886e-04 1.91e+00 8.64e-06 -5.7 6.80e+01 - 1.00e+00 5.00e-01h 2
56 5.0744550e-04 3.37e-01 1.64e-05 -5.7 2.54e+01 - 1.00e+00 1.00e+00h 1
57 5.0763613e-04 1.82e+00 1.43e-05 -5.7 5.79e+01 - 1.00e+00 1.00e+00h 1
58 5.0751252e-04 2.25e-02 1.93e-05 -5.7 6.65e+01 - 1.00e+00 1.00e+00H 1
59 5.0783220e-04 1.23e-04 1.85e-05 -5.7 6.92e+01 - 1.00e+00 1.00e+00H 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
60 5.0751252e-04 1.60e-02 1.93e-05 -5.7 6.89e+01 - 1.00e+00 1.00e+00H 1
61 5.0783207e-04 3.26e-04 1.85e-05 -5.7 6.92e+01 - 1.00e+00 1.00e+00H 1
62 5.0751786e-04 7.95e-05 1.93e-05 -5.7 6.89e+01 - 1.00e+00 1.00e+00H 1
63 5.0754291e-04 6.40e-01 8.69e-06 -5.7 6.93e+01 - 1.00e+00 5.00e-01h 2
64 5.0745167e-04 3.00e-01 1.61e-05 -5.7 2.35e+01 - 1.00e+00 1.00e+00h 1
65 5.0763437e-04 1.71e+00 1.33e-05 -5.7 5.61e+01 - 1.00e+00 1.00e+00h 1
66 5.0751762e-04 2.91e-03 1.93e-05 -5.7 6.60e+01 - 1.00e+00 1.00e+00H 1
67 5.0752104e-04 1.62e-01 1.43e-05 -5.7 6.92e+01 - 1.00e+00 2.50e-01h 3
68 5.0761295e-04 1.03e+00 2.18e-05 -5.7 4.35e+01 - 1.00e+00 1.00e+00h 1
69 5.0751251e-04 2.38e-02 1.92e-05 -5.7 5.95e+01 - 1.00e+00 1.00e+00H 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
70 5.0783210e-04 3.29e-06 1.85e-05 -5.7 6.90e+01 - 1.00e+00 1.00e+00H 1
71 5.0751786e-04 7.94e-05 1.93e-05 -5.7 6.89e+01 - 1.00e+00 1.00e+00H 1
72 5.0754292e-04 6.40e-01 8.69e-06 -5.7 6.93e+01 - 1.00e+00 5.00e-01h 2
73 5.0749299e-04 3.95e-01 1.21e-05 -5.7 2.35e+01 - 1.00e+00 5.00e-01h 2
74 5.0755524e-04 1.36e-01 3.09e-05 -5.7 1.57e+01 - 1.00e+00 1.00e+00h 1
75 5.0751308e-04 1.96e-02 1.80e-05 -5.7 3.30e+01 - 1.00e+00 1.00e+00H 1
76 5.0754492e-04 5.62e-01 8.19e-06 -5.7 6.42e+01 - 1.00e+00 5.00e-01h 2
77 5.0748301e-04 3.98e-01 1.25e-05 -5.7 2.93e+01 - 1.00e+00 5.00e-01h 2
78 5.0756455e-04 2.42e-01 3.07e-05 -5.7 2.09e+01 - 1.00e+00 1.00e+00h 1
79 5.0751710e-04 1.38e-03 1.87e-05 -5.7 3.98e+01 - 1.00e+00 1.00e+00H 1
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
80 5.0783077e-04 3.26e-06 1.84e-05 -5.7 6.70e+01 - 1.00e+00 1.00e+00H 1
81 5.0726420e-04 1.24e+01 1.85e-05 -8.6 1.34e+02 - 9.93e-01 1.00e+00h 1
82 5.0749897e-04 2.75e+00 2.53e-06 -8.6 6.93e+01 - 1.00e+00 1.00e+00h 1
83 5.0749686e-04 3.64e-04 6.40e-10 -8.6 4.69e-02 - 1.00e+00 1.00e+00h 1
84 5.0749686e-04 7.28e-11 2.51e-14 -8.6 6.86e-05 - 1.00e+00 1.00e+00h 1
Number of Iterations....: 84
(scaled) (unscaled)
Objective...............: 5.0749685809243843e-04 5.0749685809243843e-04
Dual infeasibility......: 2.5059195074646584e-14 2.5059195074646584e-14
Constraint violation....: 1.4104644499482425e-11 7.2759576141834259e-11
Complementarity.........: 2.5059035596800647e-09 2.5059035596800647e-09
Overall NLP error.......: 2.5059035596800647e-09 2.5059035596800647e-09
Number of objective function evaluations = 260
Number of objective gradient evaluations = 85
Number of equality constraint evaluations = 260
Number of inequality constraint evaluations = 0
Number of equality constraint Jacobian evaluations = 85
Number of inequality constraint Jacobian evaluations = 0
Number of Lagrangian Hessian evaluations = 84
Total CPU secs in IPOPT (w/o function evaluations) = 0.338
Total CPU secs in NLP function evaluations = 0.092
EXIT: Optimal Solution Found.
You will notice that the resulting parameter estimation problem, when using the flash unit model, will have 2952 variables and 2950 constraints. This is because the unit models in IDAES use control volume blocks which have two state blocks attached; one at the inlet and one at the outlet. Even though there are two state blocks, they still use the same parameter block i.e. m.fs.properties in our example which is where our parameters that need to be estimated exist.
Let us display the results by running the next cell.
print("The SSE at the optimal solution is %0.6f" % (obj_value))
print()
print("The values for the parameters are as follows:")
for k, v in parameters.items():
print(k, "=", v)
The SSE at the optimal solution is 0.000507
The values for the parameters are as follows:
fs.properties.tau[benzene,toluene] = -0.8987454466579063
fs.properties.tau[toluene,benzene] = 1.410449514796474
Using the data that was provided, we have estimated the binary interaction parameters in the NRTL model for a benzene-toluene mixture. Although the dataset that was provided was temperature dependent, in this example we have estimated a single value that fits best for all temperatures.
Advanced options for parmest: bootstrapping#
Pyomo’s parmest tool allows for bootstrapping where the parameter estimation is repeated over n samples with resampling from the original data set. Parameter estimation with bootstrap resampling can be used to identify confidence regions around each parameter estimate. This analysis can be slow given the increased number of model instances that need to be solved. Please refer to https://pyomo.readthedocs.io/en/stable/explanation/analysis/parmest/covariance.html#bootstrapping for more details.
For the example above, the bootstrapping can be run by uncommenting the code in the following cell:
# Uncomment the following lines
# bootstrap_theta = pest.theta_est_bootstrap(4, seed=542)
# display(bootstrap_theta)