Optimal Molecular Design of Low-Temperature Organic Fluids under

Mar 26, 2018 - Computer-aided molecular design as a mixed-integer nonlinear programming problem under uncertainty in group contribution parameters has...
2 downloads 6 Views 1MB Size
Subscriber access provided by UNIVERSITY OF TOLEDO LIBRARIES

Process Systems Engineering

Optimal Molecular Design of Low-Temperature Organic Fluids under Uncertain Conditions Oswaldo Andres-Mart#nez, and Antonio Flores-Tlacuahuac Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.8b00302 • Publication Date (Web): 26 Mar 2018 Downloaded from http://pubs.acs.org on March 27, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Optimal Molecular Design of Low-Temperature Organic Fluids under Uncertain Conditions Oswaldo Andrés-Martínez, Antonio Flores-Tlacuahuac⇤ Escuela de Ingenieria y Ciencias, Tecnologico de Monterrey, Campus Monterrey, N.L., 64849, Mexico E-mail: [email protected] Phone: +52(1) 55 4347 2804

1

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract Computer-aided molecular design as a mixed-integer nonlinear programming problem under uncertainty in group contribution parameters has been addressed. A set of new low-temperature organic compounds, for heat recovery purposes, was obtained by solving the mixed-integer nonlinear programming problem with nominal values from a previous work. Monte Carlo simulations with Latin hypercube sampling were carried out to asses the effect of uncertain group contributions on thermo-physical properties. Furthermore, a set-based robust counterpart was formulated taking into account uncertainty only in linear constraints. The results show that even small uncertainty in group contribution parameters can lead to significant variations in thermophysical properties of the compounds analyzed. Therefore, it is necessary to consider uncertainty in the Computer-aided molecular design formulation. Solutions of the robust counterpart became more conservative as the uncertainty set size increased producing organic compounds different from the nominal case.

2

ACS Paragon Plus Environment

Page 2 of 40

Page 3 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

1

Introduction

There exists a huge interest in replacing traditional energy sources (i.e. fossil fuels) because of their negative impact on the environment 1 and their finite availability. 2 Furthermore, global energy consumption continues increasing as population and economy grow. 3 Moreover, it is important to reduce the generation and emission of greenhouse gases (particularly CO2 ) so that industrial activities do not contribute remarkably to climate change. 4 There are well-known renewable and sustainable energy sources which are feasible alternatives to fossil fuels. Examples of these sources are biofuel, biomass, geothermal, hydropower, solar, wind power, and tidal power. Currently, they supply about 23.7% of the total world energy demand and many of them are still under development. 5 As a result, research in favor of the renewables sector is an ongoing concern these days. In addition to using these alternatives, recovery of energy from waste heat is an appealing way to maximize energy usage and reduce the negative environmental impact. This energy is discharged due to its low temperature and considered as useless, hence it is essential not only to address energy recovery but also to increase conversion efficiencies to reuse it effectively. Organic Rankine Cycle (ORC) is an example of a mature technology for waste heat recovery. An ORC is essentially a Rankine cycle that uses an organic working fluid instead of water. 6 It is usually applied for conversion of geothermal and biomass energy into electricity. 7,8 However, the application in waste energy recovery has gained some interest. Heat can be recovered at much lower temperature due to the lower boiling point of a properly selected organic working fluid. Lecompte et al. 9 analyze advantages and disadvantages of different architectures in which the ORC can be configured for this application. Stijepovic et al. 10 explored the relationships between working fluid properties and ORC common economic and thermodynamic performance criteria, they showed that the performance of ORC systems strongly depends on working fluid properties. When selecting the most appropriate working fluid, Quoilin et al. 7 recommend to consider a 3

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

maximum temperature of the working fluid for an optimal ORC. Other aspects to be looked at are environmental (e.g ozone depletion potential and global warming potential) and safety (e.g flammability and auto ignition) issues. 11 Therefore, it is important to find the best compound(s), either existing or new, that make energy recovery more efficient and environmental friendly. Molecular design is an approach commonly utilized to obtain new compounds or molecules with desirable properties. Computer-aided molecular design (CAMD) is a technique that has been efficiently adopted for designing optimal molecular structures for different applications. 12 Group contribution (GC) methods are used as quantitative structure-property relationships within the CAMD problem in order to estimate molecular properties from a determined set of functional groups. 13 These methods take information about the molecular structure, that is the frequency and type of a functional group appearing in the molecule. One of the GC methods employed in CAMD was proposed by Joback and Reid, 14 which is an extension of the Lydersen method. 15 They applied linear regression techniques to determine the group contributions for a set of functional groups covering a wide variety of organic compounds and assuming no interaction between groups. Accuracy and applicability of GC models have been improved over the years. 16–21 The CAMD problem is usually cast as an optimization problem, which consists of inequality and equality constrains that ensure structural feasibility and take into account thermodynamical properties and an objective function modeling a goal to be optimized. 22–25 This principle has been applied for designing optimal molecular structures to be used as organic working fluids in ORC. In this approach the CAMD problem is usually formulated as a mixed-integer nonlinear programming (MINLP) problem where binary variables determine the structure of the new molecules with desired properties. 26,27 There are some studies tackling the problem with multi-objective optimization algorithms to optimize both process performance and working fluid structure. 28,29 Optimal mixtures of 4

ACS Paragon Plus Environment

Page 4 of 40

Page 5 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

organic working fluids have been also investigated by using MINLP and multi-objective techniques. 30–32 Estimated parameters (group contributions) in GC methods are subject to uncertainty, leading to uncertain predicted properties. As a result, GC methods provide only estimates, and depending on the particular application the inaccuracy of the results can be significant and vary according to the property being calculated. 33 Hence, assessing uncertainty of both estimated parameters and predicted properties is an important key to obtain robust solutions. When uncertainty information is not reported, it is common to quantify it by regression techniques applied to experimental data and determine the type of probability distribution as well as mean and variance of the uncertainty. 34–38 Once uncertainty information is obtained or found in literature, a common strategy is propagating the uncertainty through the model with an appropriate sampling method and evaluating the output with statistical techniques. 39,40 Reed and Whiting 41 used Monte Carlo simulations with Latin Hypercube Sampling (LHS) to propagate the influence of the uncertainty in thermodynamic data on the performance of a binary distillation model. They represented the results as empirical cumulative distribution functions (CDF) and determined the individual effects calculating partial correlation coefficients (PCC). A similar approach was followed by Whiting et al. 42 applied to two distillation models. Frutiger et al. 43 carried out MC simulations with LHS to propagate the influence of the input uncertainty of fluid properties on an ORC model and calculated a 95% confidence interval for the outputs. Hajipour and Satyro 44 applied a traditional random MC sampling to propagate the uncertainty of thermo-physical properties through the Peng-Robinson equation of state. They showed that relatively small sample sizes in the order of 100 randomly distributed inputs may be adequate to provide reasonable uncertainty estimates for values calculated from complex models. In the present work we apply MC simulations with LHS to propagate uncertainty in group contribution parameters through the GC method and display the effect on the 5

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

output for different thermodynamic properties. Besides uncertainty quantification, the impact of uncertain parameters on the solution of the underlying optimization problem should be studied. CAMD problems with data uncertainty can be tackled with two major approaches from the field of optimization under uncertainty, namely stochastic optimization (based on random sampling) and robust optimization (RO). The first one is applied when uncertain inputs are modeled with known or assumed probability distributions. Santos-Rodriguez et al. 45 formulated a stochastic nonlinear optimization problem to address uncertainty in ORC variables (heat source temperature and turbine efficiency) for organic fluid mixture design. They showed that it is possible to get compositions that keep certain working fluid performance under different scenarios. Conversely, the second approach is useful when no information about probability distribution is available and uncertain parameters are allowed to take any realization in a uncertainty set. 46 Only a few works have focused on uncertainties in GC parameters for CAMD applications. Maranas 47 formulated a nonlinear stochastic version of the CAMD problem formulated in 22 utilizing probability density distributions to describe the likelihood of different realizations of group contribution parameters. The stochastic problem was transformed into a deterministic MINLP which allowed to figure out the effect that property prediction uncertainty may have in optimal molecular design. Kim and Diwekar 36 applied a stochastic optimization framework based on random sampling for solvent selection. To our best knowledge, uncertainty in group contribution parameters for CAMD problems featuring organic working fluids synthesis has not been addressed. Moreover, CAMD problems in general has not yet been addressed within a RO framework. The robust optimization concept we adopt in this work is based on the approach introduced by Soyster, 48 who reformulated a linear programming (LP) problem so as to obtain a robust counterpart (RC), whose solutions would be feasible under all possible perturbations in data. Other definitions for RO and examples of applications are presented by Beyer and Sendhoff. 49 Singh 50 extended Soyster’s work to fractional programming 6

ACS Paragon Plus Environment

Page 6 of 40

Page 7 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

problems. Ben-Tal and Nemirovski 51 reduced the conservatism of Soyster’s solutions by proposing an ellipsoidal-set-based RC. Ben-Tal and Nemirovski 52 and El-Ghaoui et al. 53,54 developed a similar approach to deal with parameter uncertainty within semidefinite and quadratic problems. They found that these robust counterparts are computationally tractable. Ben-Tal and Nemirovski 55 suggested a methodology to treat data contaminated with uncertainty either unknown-but-bounded or random symmetric in LP problems. They formulated a RC based on probability of constraint violation and infeasibility tolerance. This methodology was extended to mixed-integer linear programming (MILP) problems by Lin et al. 56 Bertsimas and Sim 57 introduced a budget parameter to control the degree of conservatism of previous formulations. Li et al. 56 developed RCs based on different shapes and size of uncertainty sets for both LP and MILP problems. RO has been applied in engineering problems mostly in planning and scheduling studies, but it has been also applied in structural design, circuit design, power control in wireless channels, antenna design, linear-quadratic control and many others. 58–60 However, RO theory for general nonlinear programming (NLP) problems is not as well established as for the linear case. Zhang 61 proposes a formulation for NLP problems that is valid in the neighborhood of the nominal values. Hale and Zhang 62 apply this RC to the design of a heat exchanger network and reactor-separator system. More recent developments in RO can be found in the review by Gabrel et al. 63 In this work, we use techniques of linear robust optimization which can be applied at constraint level to obtain a RC of the original MINLP problem. This paper is organized as follows: Section 2 presents the original CAMD technique as a MINLP problem to generate a family of new organic compounds, the strategy to carry out MC simulations with LHS for uncertainty analysis, and the RC formulation. In Section 3 we discuss the results obtained from the MINLP with nominal values, the uncertainty analysis, and solutions of the RC for different cases. Finally in Section 4 we list our main conclusions about this work. 7

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

2

Methodology

Palma-Flores et al. 26 applied CAMD to obtain a new family of organic fluids for energy recovery from low-temperature energy sources. They started from a basis set of 39 functional groups and solved several MINLP problems where binary variables determine the structure of the new compounds with desired properties. The GC method used was that of Joback and Reid. 14 We use the same formulation with a reduced basis set of functional groups to get a family of 12 organic compounds. For each of the compounds obtained with nominal values, a MC simulation with LHS was implemented in order to evaluate the impact of uncertainty in group contribution parameters on thermodynamic properties. Then, we formulate a RC of the original CAMD problem considering uncertainty only in the linear constraints. Finally, the robust solutions are compared with the base nominal cases corresponding to the four addressed objective functions considered in. 26

2.1

CAMD nominal case

Following the same approach of Palma-Flores et al., 26 a version of the CAMD problem was written with a basis set of 14 functional groups including only those appearing more frequently in the 32 molecules obtained by Palma-Flores et al. 26 These functional groups are: CH3 , CH2 , CH, C, F , Cl, OH, O, C, COO, N H2 , N H, N , and S. The optimization

8

ACS Paragon Plus Environment

Page 8 of 40

Page 9 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

problem has the following form: (1a)

(x, y)

max/min x,y

s.t 0  Ay  b

(1b)

H(y, x, a) = 0

(1c)

xL  x  xU

(1d)

y 2 {0, 1}, x 2 Rn where

is one of the four objective functions studied in. 26 Constraints (1b) define the

structural feasibility of the organic fluid. Constraints (1c) set the GC method used to estimate the physical and thermodynamic properties which are represented by the continuous decision variables x. Upper and lower limits denoted by U and L, respectively, are imposed on x. Binary variables y determine the molecular structure of the organic compound. Finally a are the GC parameters. The complete MINLP problem formulation is provided in the Supporting Information. The four objective functions are: • Maximization of enthalpy of vaporization max

1

H lv ,

=

H lv

• Maximization of the ratio between enthalpy of vaporization

(2)

H lv and liquid heat

capacity Cpl , max

2

9

=

H lv Cpl

ACS Paragon Plus Environment

(3)

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 40

• Maximization of weighted sum of properties, max where

3

= 0.5 H lv

0.1Cpl

0.4 Gf298

(4)

Gf298 is the Gibbs energy of formation at 298 K.

• Minimization of Cpl , min

4

= Cpl .

(5)

We tested the four objective functions in the original problem (Eqs.2-5) and varied the minimum (nmin ) and maximum (nmax ) number of functional groups allowed to form the compound. For example, if one resulting molecule consists of four (ten) functional groups, the same program is run with nmin = 5 (nmax = 9). Hence, it is possible to obtain molecules of different size for each objective function. A total of 12 MINLP problems (three for each objective function) were solved with nominal values of group contributions. As Tb is a key property for practical purposes, the upper bound Tb  373.15 K should be kept in mind for subsequent results discussion. We also avoided undesirable attachments like -O-OH and -COO-O- as suggested by Palma-Flores et al. 26 The optimization problems were implemented in GAMS (version 24.2.3) 64 with Dicopt as MINLP solver. 65

2.2

Uncertainty analysis

As pointed out by Maranas, 22 contribution of molecular groups for a given property is not necessarily unchanged from one molecule to another, but varies slightly around some nominal value. This behavior can be expressed assuming that group contribution parameters follow a normal distribution with mean µ and standard deviation . With the information reported by Joback 66 for each property, we take the nominal group contribution value as mean µ and, since standard deviation is not reported, we take the average 10

ACS Paragon Plus Environment

Page 11 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Table 1: Uncertainty associated with each property model. Property Tb Average per cent error 3.6

Tf 11.2

Tc 0.8

Pc 5.2

Vc 2.27

Hf 9.2

Gf 15.7

H lv 3.88

Cpg 1.4

of the absolute error associated with each GC model (Table 1) as a reference for setting . To investigate the effect of uncertainty in group contribution parameters on predicted properties, MC simulations with LHS were performed to propagate the input uncertainty through the GC method model of Joback and Reid. 14 The corresponding equations for properties in Table 1 are Eqs. S21, S22, S23, S25, S27, S28, S29, S30 and S31 from the Supporting Information. MC simulation involves three steps: (1) specifying input uncertainty (by a distribution); (2) sampling input uncertainty and (3) propagating the sample through the model. LHS displays properties between random sampling, which involves no stratification, and stratified sampling, which stratifies all the sampling space. 67 As illustrated in Figure 1, the procedure is the inverse to CAMD. We start with a fixed molecule, then for a given property we specify µ and

so that LHS generates N sam-

ples and passes them to the GC method model. The results are depicted in empirical CDFs, scatterplots and PCCs. The MC simulations were programmed in R software (version 3.4.0) 68 with the parameter space exploration (pse) package 69 with N = 200 and 50 bootstrap replicates for calculating the PCCs.

Fixed molecule

GC method model

LHS μ, σ

N samples

Figure 1: Procedure for MC simulations with LHS.

11

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

2.3

Page 12 of 40

Robust formulation of uncertain linear constraints

The RO concept adopted in this part is based on the worst case formulation introduced by Soyster. 48 Two approaches are followed, namely those in Lin et al. 56 and Li et al. 70 They obtained RCs for linear and mixed-integer linear optimization problems. Their formulations can be applied at (linear) constraint level. However, a general RO theory for the mixed-integer nonlinear case has not been developed yet. Accordingly, we treat only the linear constraints in the MINLP problem as a new way to address uncertainty in this framework. The original MINLP problem can be rewritten including uncertainty in linear equalities G of GC model (Eqs. S23, S25, S28, S29, S30 and S31 in the Supporting Information) and keeping nominal values in nonlinear equalities F as follows

min/max

(x, y)

(6a)

0  Ay  b

(6b)

F (y, x, a) = 0

(6c)

G(y, x, e a) = p

(6d)

s.t

(6e)

xL  x  xU

y 2 {0, 1}, x 2 Rn , where p are the constants appearing in each linear GC equation, e a are the values of coefficients subject to uncertainty and a the nominal values. The jth linear constraint in (6d) with uncertain coefficients has the general form

xj

nX m max X i=1 k=1

yik e ajk = pj ,

(7)

ajk is the group where xj is one of the continuous variables Tb , Tf , Vc , Gf , H f , H lv ; e 12

ACS Paragon Plus Environment

Page 13 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

contribution k subject to uncertainty, that is e ajk = ajk + b ajk ⇠jk ; ajk is the nominal value;

b ajk = ⌧j ajk represents perturbation around the nominal value given by ⌧j > 0, whose

values are the errors in Table 1; and ⇠jk are random variables distributed in the interval

[ 1, 1]. Thus, when ⇠jk =

1, e ajk = ajk

⌧j ajk ; when ⇠jk = 1, e ajk = ajk + ⌧j ajk ; and when

ajk = ajk , that is, the uncertain parameter is at its nominal value. As an example, ⇠jk = 0, e

Tb constraint in the form of (7) becomes

Tb

nX m max X i=1 k=1

yik Tebk = 198.2,

(8)

with Tebk = Tbk + Tbbk ⇠k and Tbbk = ⌧ Tbk . The goal is to formulate a RC of (6) able to find solu-

tions that remain feasible for any realization of ⇠ in a given uncertainty set. According to

Lin et al., 56 for each inequality constraint involving uncertain coefficients, an additional constraint is introduced to incorporate the uncertainty and maintain the relationships, in this case the quantitative structure-property relationships (GC method), among the relevant variables under a given infeasibility tolerance. Assuming that uncertain coefficients e ajk are randomly and symmetrically distributed around the nominal values, two approaches can be followed to derive a RC of (7).

The first approach is proposed by Lin et al. 56 which is an extension of Ben-Tal and Nemirovski 55 to mixed-integer problems. A solution (x, y) is robust if satisfies the following conditions: (i) (x, y) is feasible for the nominal problem (ii) For the jth inequality, the probability of the constraint violation is at most  > 0, with an infeasibility tolerance > 0, where  is a reliability level

13

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 40

Since the equality constraint in (7) is equivalent to the pair of inequalities 71 nX m max X

xj xj +

i=1 k=1 nX m max X i=1 k=1

yik e ajk  pj yik e ajk 

(9a) (9b)

pj ,

the RCs of (9a) and (9b) based on Lin et al. 56 (Theorem 2) are respectively nX m max X

xj

i=1 k=1

xj +

nX m max X i=1 k=1

v unmax m uX X yik ajk + ⌦j t yik b a2jk  pj +

j

i=1 k=1

v unmax m uX X yik ajk + ⌦j t yik b a2jk  i=1 k=1

pj +

max{1, |pj |}

j

max{1, |pj |},

(10a)

(10b)

where ⌦j is a positive parameter, and (x, y) is a robust solution that satisfies (i) and (ii) with j = exp

⌦2j /2 .

ajk in terms of ⇠jk in (9a)-(9b) The second approach comes from Li et al. 70 Replacing e

and grouping, a worst-case formulation is obtained:

xj xj +

nX m max X i=1 k=1 nX m max X i=1 k=1

nX m max X

yik ajk + max{ ⇠j 2⌅

i=1 k=1 nX m max X

yik ajk + max{ ⇠j 2⌅

i=1 k=1

yik ⇠jk b ajk }  pj yik ⇠jk b ajk } 

(11a) pj .

(11b)

There are many options to choose a set ⌅ where ⇠jk can take any realization. However, as we have assumed e ajk are randomly and symmetrically distributed around the nominal

values, the appropriate set to be chosen is an ellipsoidal one, since a box would be too pessimistic. 52,72 An ellipsoidal set is defined with 2-norm as follows:

⌅ = {⇠|k⇠k2  where

}

(12)

is the adjustable parameter controlling the size of the uncertainty set. Thus, the 14

ACS Paragon Plus Environment

Page 15 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

RCs of (11a)-(11b) induced by the set ⌅ are nX m max X

xj

i=1 k=1

xj +

nX m max X i=1 k=1

v unmax m uX X yik ajk + j t yik b a2jk  pj + v unmax m uX X yik ajk + j t yik b a2jk 

where an infeasibility tolerance

j

i=1 k=1

pj +

j

i=1 k=1

j

(13a)

max{1, |pj |}

(13b)

max{1, |pj |},

has been added. Note that when ⌦j =

j,

both for-

mulations (10) and (13) are equivalent. This is due to the fact that in the original problem uncertain coefficients take place only in the binary variables, otherwise the RCs would have been slightly different. In the following, we take (10) and (13) as the same formulation with parameter ⌦j . Consequently, the complete RC of the CAMD problem (6) is (x, y)

(14a)

0  Ay  b

(14b)

F (y, x, a) = 0

(14c)

G(y, x, a) = p

(14d)

max/min x,y

s.t

xj

v unmax m nX m max X uX X yik ajk + ⌦j t yik b a2jk  pj + i=1 k=1

i=1 k=1

v unmax m nX m max X uX X xj + yik ajk + ⌦j t yik b a2jk  i=1 k=1

i=1 k=1

xL  x  xU

j

max{1, |pj |},

pj +

j

(14e)

8j

max{1, |pj |},

8j

(14f) (14g)

y 2 {0, 1}, x 2 Rn . The inclusion of constraints (14e)-(14f) means incorporating the worst case values b ajk 15

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 40

into the MINLP problem. Hence, these constraints are the most difficult to maintain, so the parameter

allows a tolerable violation. As a consequence, binary variables y

will now search for an optimal combination of functional groups whose contributions satisfy the GC models (Eqs. (14c)-(14d)) while taking into account uncertainty through constraints (14e)-(14f). As an example, consider the GC method constraint for Tb (8). The additional constraints are

Tb

nX m max X i=1 k=1

Tb +

nX m max X i=1 k=1

v unmax m uX X 2 yik Tbk + ⌦Tb t yik Tbbk  198.2 (1 +

Tb )

(15)

1)

(16)

i=1 k=1

v unmax m uX X 2 yik Tbk + ⌦Tb t yik Tbbk  198.2 (

Tb

i=1 k=1

The RC of the CAMD problem was solved for each objective function with  = 10% (⌦j = 2.1459),

j

= 5% for all j except for

Gf , V c y

H f where ⌦

Gf

= ⌦V c = ⌦

Hf

=

0.1⌦j . Solutions of the nominal case were taken as initial values for the respective uncertain optimization cases keeping the same bounds. In this part, we are not interested in varying the minimum and maximum number of functional groups allowed to form the compound. Therefore, we take only the base nominal case for each objective function as a reference for the RC. In addition, for each objective function the problem (14) was solved for different values of ⌦j and

= 5%, and finally a sensitivity analysis with respect to

j

was done so as to observe its effect on the robust solution.

3 3.1

Results and discussion CAMD nominal case

Table 2 lists the compounds obtained with nominal values in group contribution parameters for each run. The objective values along the boiling temperature are also shown. Solutions 1, 4, 7 and 10 are the base cases, i.e. the optimal solutions found by the MINLP 16

ACS Paragon Plus Environment

Page 17 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

problem without any additional restrictions. Solutions 2 and 3 are the result of varying the minimum and maximum number of functional groups allowed to be part of the molecular structure of base case 1. Solutions 5 and 6 are the result of such a variation in base case 4. Similarly, solutions 11 and 12 come from the same variation in base case 10. When limitations on the number of functional groups are imposed, the value of the objective function decreases for

1,

2

and

3

but increases for

4

since it is a minimiza-

tion problem. In compounds 1, 2 and 3 one attachment N (OH) is part of the molecular structure. In compounds 4, 5 and 6, functional group Cl appears twice. Functional group COO appears only in compounds associated with objective function

3

and takes place

just once in them. Functional group F is present in 9 out of 12 compounds. From Table 2, only compounds 1 (CH3 -O-N(OH)-CH3 ) and 4 (Cl-O-S-Cl) are reported by Palma-Flores et al. 26 Compounds 10 and 11 are found in literature as chlorine peroxide and 1,1-dichloro1-fluoroethane respectively. The highest Tb belongs to compound 1 and the lowest to compound 11. All the compounds in Table 2 are candidates to be utilized as working fluids for heat recovery at low temperature, since their boiling temperatures are below that of water. Table 2: Compounds obtained with nominal values Objective Function and optimal value

ID

Compound

1 2 3 4 5 6 7 8 9 10 11 12

CH3 O N (OH) CH3 F O O N (OH) CH3 F N (OH) CH2 O O F Cl O S Cl F S N (Cl) Cl F C(Cl, F ) S Cl F COO C(O F )3 F COO N H C(F )2 O F F COO C(F )2 S F Cl O O Cl CH3 C(Cl)2 F F S C(F )2 Cl

17

1 1 1 2 2 2 3 3 3 4 4 4

ACS Paragon Plus Environment

= 41.178 = 40.545 = 39.728 = 0.2765 = 0.2313 = 0.2010 = 565.02 = 446.77 = 429.44 = 115.49 = 133.70 = 142.04

Tb (K) 372.40 371.20 370.48 365.66 354.95 361.43 364.69 370.02 366.21 319.30 316.26 323.27

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

3.2

Uncertainty analysis

Uncertainty analysis was performed for all compounds obtained with the nominal case of the CAMD problem (Table 2) and all the properties in Table 1. The empirical CDF (Figures 2(a) and 3(a)) is basically a step function that jumps up 1/N at each data point. It estimates the underlying CDF generated by the points in the sample. As an example, consider the compound CH3 -O-N(OH)-CH3 . Figure 2 shows the results of MC simulations for Tb . It can be observed in the empirical CDF (Figure 2(a)) that for the given variation in group contributions, Tb varies from 362 K to 381 K, which can represent an important uncertainty for CAMD purposes as Tb is required to be less or equal than 373.15 K. The histogram drawn in the same figure provides a sight of this variation. When averaging 200 model outputs, the average Tb is 372.4 K, which is the Tb calculated with nominal values of GC parameters. The PCC (Figures 2(b) and 3(b)) measures the strength of the linear association between the output and each GC after removing the effect of the others. This relation can be visualized in the scatterplots (Figures 2(c) and 3(c)) that show the distribution of values returned by the GC model and how sensible these responses are with respect to variations in the values (x-axes) of each GC parameter. As shown in Figure 2(b) OH group has the strongest effect on Tb for this particular molecular architecture when removing the effect of the other groups, whereas N group has the weakest effect. The points in the scatterplot (Figure 2(c)) Tb versus OH contribution lie near the straight line. Uncertainty analysis for Pc is depicted in Figure 3. In this case the corresponding GC method is nonlinear. The empirical CDF along with the histogram in Figure 3(a) suggests little variation in Pc respect to uncertainty in group contribution parameters, namely 54.20 bar to 57.82 bar. PCCs in Figure 3(b) show that N and OH groups have a strong linear effect on Pc being OH contribution slightly stronger. In this case, O group has the weakest effect. This can be viewed in the scatterplots of Figure 3(c). The average Pc for the compound under analysis is 55.95 bar, same as calculated from nominal values. This graphical analysis for all the compounds and thermo-physical properties, as well as 18

ACS Paragon Plus Environment

Page 18 of 40

Page 19 of 40 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

the R-file are available in the Supporting Information. The resulting Tb of the MC simulations for the 12 compounds are depicted as boxplots in Figure 4(a). Each box is drawn by the first quartile, median (thick vertical line) and third quartile. The whiskers represent the values which are outside the box, but within 1.5 times the interquartile range. Extreme values being outside the box are represented as small circles. Variation in Tb is similar for all the compounds. Some boxplots present extreme values on both sides. The boxplot of compound 1 has only one point outside. Compounds 8, 9, 10 and 12 have no extreme values. Figure 4(b) contains the boxplots for Pc . In this case, the width is less than those of Tb , even though the error for Pc is higher. This can be due to the non-linearity of its model. Compound 1 has two extreme values, each at one side of the box. Compounds 4 and 10 have no points outside. The previous analysis provides information about the behavior of GC method under uncertain GC parameters. Some properties can be very sensitive to input uncertainty. Therefore, it is important to take into account uncertainty within the CAMD problem since the effect of uncertain parameters propagates through the GC method and influences the molecule properties. As a result, molecular structure is likely to be different for each realization of uncertain parameters.

3.3

Robust counterpart

Results of the RC for each objective function with  = 10% and

j

= 5% are shown in

Table 3. The respective base nominal cases 1, 4, 7 and 10 from Table 2 are also shown for comparison. It can be noticed that for

1,

2

and

3

the RC leads to an objective value

less than the respective nominal case, with a more notable decrease for

3.

For

4,

the

objective value of the RC is greater because it is a minimization problem. Every robust version finds a well known compound with a smaller boiling temperature, and each one is different from its nominal version. Every compound resulting from the RC is feasible for any realization of the uncertainty in the set ⌅ with a given value of ⌦j . The first RC 19

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

0.5 Tb (K) 0.0

0.6

−0.5

0.4

−1.0

0.2

Proportion