Toward Optimum Working Fluid Mixtures for ... - ACS Publications

Jul 17, 2013 - evaluated using an ORC process model in the course of molecular mixture ... performance, few ORC plants utilizing working fluid mixture...
0 downloads 0 Views 2MB Size
Subscriber access provided by UNIV OF NEBRASKA - LINCOLN

Article

Towards Optimum Working Fluid Mixtures for Organic Rankine Cycles using Molecular Design and Sensitivity Analysis Athanasios Papadopoulos, Mirko Zoran Stijepovic, Patrick Linke, Panos Seferlis, and Spyros S Voutetakis Ind. Eng. Chem. Res., Just Accepted Manuscript • Publication Date (Web): 17 Jul 2013 Downloaded from http://pubs.acs.org on July 17, 2013

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 free 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 accessible to all readers and 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.

Industrial & Engineering Chemistry Research 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 61

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

Towards Optimum Working Fluid Mixtures for Organic Rankine Cycles using Molecular Design and Sensitivity Analysis Athanasios I. Papadopoulosa,* , Mirko Stijepovicb, Patrick Linkeb, Panos Seferlisc, Spyros Voutetakisa a

Chemical Process and Energy Resources Institute (CPERI), Centre for Research and Technology Hellas (CERTH), P.O. Box 60361, Thermi-Thessaloniki 57001, Greece.

b

Chemical Engineering Department, Texas A&M University at Qatar, P.O. Box 23874, Doha, Qatar c

Department of Mechanical Engineering, Aristotle University of Thessaloniki, P.O. Box 484, 54124 Thessaloniki, Greece

*Email: [email protected], Tel:+30 2310 498363, Fax:+30 2310 498360.

ACS Paragon Plus Environment

1

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 2 of 61

ABSTRACT This work presents a Computer-Aided Molecular Design (CAMD) method for the synthesis and selection of binary working fluid mixtures used in Organic Rankine Cycles (ORC). The method consists of two stages, initially seeking optimum mixture performance targets by designing molecules acting as the first component of the binaries. The identified targets are subsequently approached by designing the required matching molecules and selecting the optimum mixture concentration. A multi-objective formulation of the CAMD-optimization problem enables the identification of numerous mixture candidates, evaluated using an ORC process model in the course of molecular mixture design. A nonlinear sensitivity analysis method is employed to address model-related uncertainties in the mixture selection procedure. The proposed approach remains generic and independent of the considered mixture design application. Mixtures of high performance are identified simultaneously with their sensitivity characteristics regardless of the employed property prediction method.

KEYWORDS

Organic Rankine cycle, working fluid mixtures, molecular design, sensitivity analysis, multi-objective optimization

ACS Paragon Plus Environment

2

Page 3 of 61

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

Introduction The utilization of low-grade heat sources for power generation has received significant attention in the past decade in view of increasing concerns over climate change and high energy prices. Organic Rankine cycle (ORC) systems constitute promising technologies to exploit widely available resources of low thermal content. Their operation is based on the vaporization of a working fluid to drive a turbine, hence there are two major characteristics affecting the thermodynamic and economic cycle performance for a given heat source; the type of the employed working fluid and the ORC system features. ORCs commonly utilize single working fluids to facilitate efficient heat extraction. Existing reports indicate that their replacement by binary or multi-component organic mixtures in conventional ORCs1 or regenerative systems2 may significantly increase their efficiency, while the levelized electricity costs and cooling expenditures decrease by 24%3 and 40%4 for geothermal and waste heat recovery applications. With the potential to also improve the environmental and safety system performance, few ORC plants utilizing working fluid mixtures are also in operation5. Mixtures enable higher overall system performance compared to single fluids for various reasons. Their evaporation temperature is adaptively adjusted to the profile of the heat source, as opposed to the flat evaporation profile of pure fluids. This results in a condition coined “Temperature glide”1 which occurs due to the continuous change in the liquid composition of the mixture at boiling conditions, generating a variable mixture boiling temperature. Temperature glide avoids the appearance of a pinch as in the case of pure fluids and significantly reduces exergy losses. Furthermore, mixtures introduce design flexibility as they combine favorable environmental and safety properties with increased economic and operating performance through regulation of the mixture composition. The consideration of different working fluid mixtures for ORCs has been the subject of an increasing number of publications1-15. Previously published work1 has investigated different types of multicomponent mixtures comprising hydrocarbons, hydrofluorocarbons or siloxanes together with important mixture performance measures and constraints that need to be considered for their evaluation. Hydrocarbon mixtures were also proposed considering regenerative preheating ORC schemes2, 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 4 of 61

equipment sizing3 and zeotropic fluids6 for efficient exploitation of moderate temperature geothermal resources. Mixtures of siloxanes or hydrocarbons have been considered4 to recover wasted heat from molten carbonate fuel cells using an ORC. Halocarbon mixtures have been investigated7 for power generation using geothermal heat, indicating significant ORC performance gains compared to pure fluids. A binary mixture of fluorocarbons has been investigated8 at different concentrations employed in an ORC system for power generation from solar energy. In a similar context, the utilization of hydrocarbon and fluorocarbon mixtures has been investigated9 at different temperature heat sources in ORCs, evaluating the resulting performance gains using ORC operating parameters like inlet/outlet volume ratio, mass flow, enthalpy difference of expansion etc. Different combinations of binary and tertiary mixtures have also been evaluated10 including alkanes, fluorinated alkanes and siloxanes aiming to find their optimum concentration. An investigation of organic, ammonia-water and alcohol-water mixtures was performed using an optimization method to identify their optimum concentration in ORC and Kalina cycle systems11. Mixtures of ammonia- water and CO2- water were also considered in two new ORC configurations, namely the ORC with liquid-flooded expansion and the ORC with solution circuit12, with the mixtures employed in the second configuration only. An ammonia-water mixture was also considered in the context of a Kalina cycle and its performance was compared with an ORC using pure ammonia or R134a13. Nineteen binary working fluid mixtures were also considered as an alternative to ammonia-water, resulting in the conclusion that the highest performers were propane and propylene-based mixtures. Binary and tertiary polysiloxane mixtures are considered in a different work14 for ORCs recovering heat from cogeneration plants fed with wood residuals. A zeotropic mixture of R227ea/R245fa is analysed in a subcritical ORC employed for exploitation of geothermal resources15 To enable optimum performance in ORCs, existing research efforts investigate combinations among common pure working fluids aiming to determine their optimum concentration within the mixture. Although useful, this approach is limited by the consideration of fluids drawn from empirically compiled repositories, often containing well-known, conventional options. The included fluids represent a very small set in view of the vast number of molecules that can be considered as components of a ACS Paragon Plus Environment

4

Page 5 of 61

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

potential mixture, hindering the opportunities for significant performance gains. Such challenges can be addressed using a Computer Aided Molecular Design (CAMD) approach, a technology first implemented by the authors for the design of pure ORC working fluids16 considering different properties17 illustrated in geothermal power generation18 and heat/power cogeneration19 applications. Recently, a CAMD-based approach for pure ORC working fluid design has also been proposed20 by other authors. CAMD is based on the experimentally measured contributions of chemical functional groups to a particular molecule (group contribution methods), aiming to enable the calculation of pure component (e.g. boiling point) and complex mixture properties (e.g. thermodynamic equilibria etc.). Utilizing molecular characteristics as performance measures within an optimization-based formulation it supports the synthesis of novel molecular structures or conventional but previously overlooked molecules with optimum features in terms of both properties and process performance. Compared to the wide use of refrigerants or hydrocarbons in existing ORC systems, the single fluids designed using CAMD16 indicate increased operating and economic ORC performance, improvements in safety and environmental properties, intense similarities with commercially patented working fluids as well as a previously unexplored molecules containing a broad range of combinations among amine-, fluorine-, ether- and other groups. The above benefits constitute a significant motivation for the development and implementation of a CAMD-based method for the design of ORC working fluid mixtures. However, considerable challenges are associated with such an effort due to the need for simultaneous determination of the optimum number of mixture components, composition (i.e. chemical structure of all participating working fluids) and concentration (i.e. amount of each fluid in the mixture). The resulting combinatorial complexity requires intense computational effort, which is further enhanced by the need for ORC models in the course of CAMD optimization to predict complex mixture behaviours (e.g. enthalpy changes, vapourliquid equilibrium etc.). CAMD-based mixture design has yet to be considered in ORC research, while there is a number of reports regarding applications in the design of solvents21-26, refrigerants27-29, polymers30-32 and gasoline blends33, 34: ACS Paragon Plus Environment

5

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 6 of 61

Previously published work21 has addressed the CAMD-based design of solvent blends for environmental impact minimization considering numerous important mixture property constraints (e.g. binary azeotropic points, solvent compatibility and so forth). The method is applied for the simultaneous determination of composition and concentration of binary mixtures in absorption processes. A decomposition-based approach has also been proposed22 to address the computational complexities involved in the design of solvent blends comprising several stages that represent the gradual evaluation of constraints such as chemical feasibility, pure component properties and mixture properties and so forth. The presented applications consider a binary solvent blend where water is pre-specified as the first solvent, while the co-solvent is designed. A similar method23 has also been proposed by the same authors. However a more efficient approach, namely interval-based optimization, is used to address the non-convexities of the thermodynamic prediction models. In a different work24,

25

the design of solvent mixtures is addressed starting from a fixed number of

components and progresses to optimize both the number of components in the mixture and the components themselves. The reported application addresses the design of a binary mixture to optimize the solubility of ibuprofen, where the mixture components are finally selected from a list of pre-specified solvents acting as decision parameters. While all the above works utilized deterministic optimization methods, a genetic algorithm has also been employed26 to address the design of multi-component mixtures through a decomposed approach applied in the design of binary solvent blends for extractive distillation processes. •

In the area of refrigerants, the design of environmentally benign refrigerant mixtures is considered27 using a database of molecules. The mixture candidates are used as integer decision parameters in an optimization problem formulation that considers binary mixtures of concentration that is additionally regarded as a continuous decision parameter. A similar approach is applied28 on a double-evaporator refrigeration cycle to obtain mixtures that maximize cooling. A different work29 focuses on addressing the combinatorial complexity in view of different concentration requirements

ACS Paragon Plus Environment

6

Page 7 of 61

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

in multi-component refrigerant mixtures considering a mixture of four refrigerants commonly used in the liquid natural gas. •

In the case of polymer blend design the development of binary polymer mixtures has been addressed30 using the structure of both components and their composition as decision parameters. In a different approach31 a property clustering method is employed to determine the optimum number of chemical constituents delivering a desired property. A similar method applied on the identification of optimum methanol blends considering three components was also recently proposed32.



Finally, a decomposition-based approach has been proposed for the design of gasoline blends33, 34. The method considers a large number of potential pure component combinations which are gradually eliminated using property performance criteria. This results in three bio-derived chemicals forming mixtures with gasoline components. It appears that the simultaneous design of all components in binary or higher mixtures is addressed in

few existing CAMD-based contributions. On the other hand, uncertainty in the context of mixture design and selection is an additional important issue which deserves attention. Potential over- or underestimation of mixture performance during CAMD due to the employed group contribution (GC), thermodynamic and/or ORC process models may impact both on the type of the designed molecular structures and their operating features. The actual operating behavior of the designed mixtures may differ from the predictions resulting from the models employed during mixture design. While several mixtures might be less sensitive to differences between predicted performance and practical implementation, others are expected to significantly deviate from their calculated performance, eventually failing to meet the desired specifications. The issue of uncertainty has been addressed in CAMD-based pure solvent design using Stochastic Annealing35 and Genetic Algorithms36, in polymer design37 and in the design of pure solvents for optimum reaction rate constants38. An uncertainty quantification method for property predictions through GC was also recently proposed39. This work addresses for the first time the CAMD-based design of binary working fluid mixtures for ORC systems. This is approached through a new method that enables the simultaneous determination ACS Paragon Plus Environment

7

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 8 of 61

(synthesis) of the working fluid structures and their concentration within the mixture, while also addressing the sensitivity involved in the selection of the optimum mixtures. The design problem is decomposed in two stages assembled around the use of a CAMD technology employing a multiobjective optimization (M.O.O.) formulation to determine mixtures with a broad range of optimum operating and safety features. The impact of uncertainty is also addressed through a systemic, non-linear sensitivity analysis approach. The proposed approach aims to a) identify the model parameters that mainly affect the ORC performance of the designed mixtures, b) quantify their sensitivity in variations imposed on these parameters through an inclusive index and c) incorporate the procedure in the selection decision-making after the implementation of CAMD. Diverse potential realizations of the predictions obtained for important properties through different and widely utilized GC methods are also accounted for during CAMD in the presented case study. This serves to impose the effects of uncertainty in mixture design to a certain extent, but the systemic consideration of uncertainty during CAMD-based mixture design is beyond the scope of this work. Overview of multi-objective CAMD for pure working fluids The design of molecules using CAMD is based on the systematic combination of molecular (functional) groups with the aim to synthesize a molecule of particular chemical structure and physicochemical properties. Such properties are calculated using GC methods, which are based on databases containing a pre-registered contribution of each group in the molecular structure containing this group. The CAMD method employed for the design of single ORC working fluids follows a systemic molecular representation approach40 which is combined with a multi-objective optimization problem formulation41 and a statistical analysis approach for integration with process design42. Additional details on the chemical principles and the operations performed to identify different molecules based on particular chemical features are reported in a previous publication40. This work only briefly describes the formal mathematical representation of molecules employed within the CAMD approach40 which is further used to develop the proposed mixture design method. Molecules are described as a set of functional groups allowed to link together based on connectivity constraints ACS Paragon Plus Environment

8

Page 9 of 61

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

imposing chemical principles to establish feasible structures. Formally a molecular structure can be described as follows: Let G = { gl }l =1 , k ∈ N be the set of available functional groups. Also let k

mk = [ g ]l =1 ∈ G k be a group vector and let Ak = diag (nl )lk=1 , nl ∈ N be a composition matrix. Then a k

molecular vector" M k " can be defined as follows40: " M k " = mk ⋅ Ak

(1)

Using equation (1) the overall attainable space " M " can be defined as: "M " =

U 1≤ k ≤ kmax

{" M

k

":" M k " = mk ⋅ Ak , mk ∈ G k , Ak = diag (nl )lk=1 , nl ∈ N , l = 1, 2,...k

}

(2)

The quotation marks are used to indicate that the considered vectors may contain both feasible or infeasible molecular structures and this is further elaborated below. The molecular structures represented through " M k " are optimized41 against a desired performance measure using an optimization algorithm such as Simulated Annealing (SA). The performance measure may result from calculations involving a detailed model of the process system hosting the designed molecules or from simpler calculations of properties capturing major operating drives of the process system. In the context of multi-objective optimization41, the employed performance measure involves a set of indices (i.e. objective functions) representing design targets. In mathematical terms, the multi-objective optimization CAMD problem is formulated as follows: Let a vector of Nof objective functions "F" =  F j 

N of j =1

forming an objective space "F" ⊆ R

N of

. Also let

" D " = [ d i ]i =d1 be a vector of Nd design variables forming a decision space " D " ⊆ R N d and let N

h ( x, d ) , q ( x, d )∀d ∈ " D ", x ∈ X be vectors of equality and inequality constraints, with X being a vector of state variables. Also let vector D ⊆ " D " including elements dj that satisfy all constraints of the form

h ( x, d ) , q( x, d ) hence representing the set of decision parameters that result in feasible molecular structures. For each d ∈ D there exists a point in the objective space "F" corresponding to mapping

ACS Paragon Plus Environment

9

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

R Nd → R

N of

Page 10 of 61

such that the space containing feasible objectives F is the image of D 43. The multi-

objective optimization problem takes the following form: min d

s.t.

F1 ( x, d ) ,..., FNof ( x, d ) h ( x, d ) = 0 q ( x, d ) ≤ 0

(3)

d L ≤ d ≤ dU x L ≤ x ≤ xU

∀d ∈ D, x ∈ X

Generally, D may represent decision variables associated with the ORC model considered in this work. Focusing in the case of molecular design it is assumed that D ≡ M , hence vector M is considered to include all feasible structures with M ⊆ " M " . Elements of M practically represent elements of vector mk and matrix Ak based on (1). In other words, changes in Mk during optimization result from replacing

the groups contained in mk by alternative ones, reducing or expanding the vector dimension k as well as altering their frequency of appearance in Ak . As gl of mk is replaced by a different group in the course of optimization, M k will represent a chemically feasible structure in M provided that equality and inequality constraints of the form h(x,d) and q(x,d) are satisfied, while their violation will result in an infeasible molecule. Indices L and U represent upper and lower bounds utilized for all the variables. Following equation (3) a feasible point d ∈ D opt ≡ M opt is called a Pareto optimum or non-dominated solution iff there exists no other point d ∗ ∈ D satisfying the following condition43: F(d ∗ ) ≤ F(d ) ∧ ∃j ∈ {1,...N of } : F j (d ∗ ) < Fj (d )

(4)

In the context of (4) the generated designs are evaluated in terms of optimality based on comparison of the objective function values representing one solution in the non-dominated set with the objective function values of the other solutions contained in the set. Additional details can be found in a previous publication41. To enable the simultaneous assessment of the desired set of objective functions through the algorithmic operations of SA, the objective function that is optimized takes the form of an aggregate objective function as follows: ACS Paragon Plus Environment

10

Page 11 of 61

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

F

Industrial & Engineering Chemistry Research

agr

N of

N in

j =1

m =1

( x, d ) = ∑ w j F j ( x, d ) + ∑ Penm qm ( x, d )

(5)

where w j is a set of weights imposed in each objective function F j , varying randomly throughout the performed optimization to enable a search of all the potential combinations41. Penm represents penalty weights that might be imposed to the inequality constraints, while N in and N of represent the total considered inequality constraints and objective functions. The continuous assessment of equation (5) through the SA algorithmic operations results in the iterative generation of new solutions. Since usually there is no single solution that optimizes (i.e. minimizes or maximizes) all the objective functions simultaneously, a set of solutions is required to be identified where at least one of the objectives is better than the others based on (4). This is the non-dominated set of solutions recorded during the search in the form of an archive which is updated based on pre-specified archival rules41.

Proposed method for mixture design The previously presented CAMD approach for the design of pure working fluids incorporates useful features that can be readily maintained in a CAMD approach targeting the design of binary mixtures. Similarly to the design of pure fluids, multiple indices need to be considered to evaluate the performance of mixtures identified through CAMD, hence the multi-objective optimization formulation is clearly a useful feature. The calculation of the performance indices may also result either from a detailed process model used in the course of optimization or from simpler property calculations capturing the influence of primary or more complex properties such as vapor-liquid equilibrium. In this respect, the existing CAMD approach for the design of pure fluids may be enhanced in a “straight forward” way by considering an additional vector of functional groups for the second component, together with the corresponding composition matrix and the mixture concentration as decision parameters to enable mixture design. Although this would result in useful solutions there are shortcomings associated with such an approach, especially when considered through a multi-objective optimization perspective: a)

The set of binary mixtures obtained as a result of multi-objective optimization would be Pareto

optimum however the chemical feasibility constraints imposed in both molecular structures in the 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

Page 12 of 61

“straight forward” approach introduces a limit in the mixture performance. Compared to the newly proposed method this may result in loss from the Pareto front of individual molecules as components of mixtures enabling high performance. A detailed analysis will be presented in the following sections as evidence for this argument. b)

The extension of the “straight forward” approach to the design of tertiary mixtures would follow

the same pattern as before, considering as decision parameters a third set of functional groups and a composition matrix with a second mixture concentration parameter. While this approach remains computationally manageable for binary mixtures, it is expected to have a detrimental effect on the computations for a tertiary or higher mixture. Although this is a likely case, the design of tertiary or higher mixtures is not addressed in this work. c)

The obtained designs in either binary or higher mixtures would have to go through a post-

processing evaluation stage to address the issue of uncertainty involved in the models employed to design the mixtures. For example, GC methods provide property predictions which may deviate from the experimentally observed values of mixtures or pure fluids, affecting also the predictions of the employed thermodynamic or ORC models. A second iteration would involve the identification of such mixtures with the aim to restore accuracy in properties or models using experimental data and to ultimately identify mixtures that are both optimum and practically useful. In other words, even if a single “straight forward” design stage results in optimum working fluid mixtures, they will not necessarily be practically applicable unless their sensitivity to variability observed in prediction methods is evaluated in an additional iteration. The merits of multi-objective optimization can be better exploited by decomposing the mixture design problem into two stages addressing the synthesis of optimum working fluids and their optimum concentration in the mixture, followed by investigation of the sensitivity of the proposed solutions to potential model variability. The following sections provide a detailed analysis of the proposed developments.

ACS Paragon Plus Environment

12

Page 13 of 61

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

Design approach Stage 1: Mixture screening For the case of a binary mixture examined in this work, the aim of stage 1 is to determine optimum molecular structures for the first component. This is approached by searching for chemically feasible molecular structures only for one of the two components, while emulating the mixture behaviour of the remaining component within a much wider structural design space by removing the chemical feasibility constraints. This serves to evaluate the performance of the first component within a broad range of potential concentrations and interactions with the second component prior to resulting in an inclusive set of molecular structures and properties for the two components. The identification of multiple optimum mixture candidates is accomplished through a multi-objective formulation of the CAMD-optimization problem, treating multiple performance measures simultaneously and resulting in a comprehensive Pareto front revealing useful structural and property trade-offs among mixture components. Figure 1 illustrates the main algorithmic steps of the proposed approach. In mathematical terms, two molecular vectors, namely M a and " M b " now represent the space of attainable mixtures " M " = M a ∪ " M b " together with a concentration vector Y = { y j } with y j ∈ [0,1] ⊆ R . By defining the concentration as a fraction of mass it holds that

∑y

j

= 1 for a mixture of Nc components, hence for a

j∈Nc

binary mixture only y1 ∈ [0,1] ⊆ R needs to be considered. Note that the key feature in the proposed approach is that molecules considered in vector M a satisfy the equality and inequality constraints of equation (3) while molecules considered in " M kb " are allowed to violate them. In other words, equation (3) is now defined for every d ∈ " D " ≡ " M " while the equality and inequality constraints are defined only for every d ∈ D ⊆ " D " . This allows the performance evaluation of more structures " M kb " generated during optimization within the same design space also used for M ka . The employed multiobjective optimization method is also key to the proposed approach as it results in a comprehensive set of mixtures " M opt " = M a ,opt ∪ " M b ,opt " containing feasible molecules in M a ,opt and feasible or ACS Paragon Plus Environment

13

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 61

infeasible molecular structures in " M b ,opt " . This is in contrast to a single optimum mixture that would be obtained in the case of single-objective optimization. Infeasible structures indicate dominating trends in terms of functional groups for the second “molecule”, while the resulting mixture properties can be used to set performance limits for the subsequent design stage. The feasible molecules that would otherwise be obtained as optimum solutions if feasibility constraints were applied in " M kb " are also evaluated since M b ⊆ " M b " , hence they are part of the available design space. In other words, although the conditions of equation (3) are allowed to be violated for " M kb " , they may also be satisfied even though this is not imposed in the optimization problem. As a result no solution is excluded from the search. Occasionally, feasible structures are outperformed by infeasible ones but since the optimization search is random the identification of feasible structures in " M b ,opt " as optimum solutions is also possible. A potential Pareto front resulting from implementation of stage 1 is illustrated conceptually in Figure 2a. It is assumed that there are two objectives (f1 and f2) that need to be maximized hence the figure

{

indicates the optimum molecules contained in M a ,opt = M ka,i

{{

}

the binary mixture contained in " M b ,opt " = " M kb,i , j "

N m1

i =1

}

N m2 j =1

}

N m1

i =1

with their corresponding structures in

, where N m1 = M a ,opt and N m 2 = " M b,opt "

are the cardinalities of the two sets. As a result, each dot point represents mixtures of the form

{

}

M a ,opt ∪ " M b ,opt " = ( M ka,1 ," M kb,1,1 "),..., ( M ka, N m1 ," M kb, Nm1 , N m 2 ") . This is a version of a Pareto front where the symbols in quotation represent only infeasible molecular structures. Notice that several feasible molecules may appear in " M b ,opt " (e.g. M kb,2,1 , M kb,2,2 etc.) while infeasible structures are also expected to appear (e.g. " M kb,1,1 "," M kb,1,2 " ). Clearly, molecule M ka,1 is part of the Pareto front in mixtures with the infeasible structures " M kb,1,1 "," M kb,1,2 " , while M ka,2 is part of the Pareto front in mixtures with feasible molecules M kb,2,1 , M kb,2,2 .

ACS Paragon Plus Environment

14

Page 15 of 61

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

Assuming that all infeasible structures of the form " M kb,i , j " (e.g. " M kb,1,1 "," M kb,1,2 " ) were to be completely removed from the Pareto front of Figure 2a, this would require the simultaneous removal of the corresponding feasible structures of the form M ka,i (e.g. M ka,1 ) since mixtures of the form

( M ka,i ," M kb,i , j ") could not be defined in such a case. Hence feasible structures of the form M ka,i (e.g. M ka,1 ) would be completely missed. On the other hand, mixtures of the form ( M ka,i , M kb,i , j ) (e.g. ( M ka,2 , M kb,2,1 ) etc.) would still remain in the Pareto front. It would also hold that no other mixture consisting of feasible structures exists that satisfies condition (4) (note that since this is a maximization problem the inequality signs of (4) must be changed appropriately). Hence there would be no other mixture covering the performance area of the missed mixtures (e.g. in the case of M ka,1 very high values for f1 at low values for f2). This is simply because if such a mixture existed it would have already been captured as a result of M b ⊆ " M b " , always assuming that the employed optimization method guarantees every time the identification of the true Pareto front. Under this assumption, if the mixture design problem was solved by implementing the chemical feasibility constraints of equation (3) in both M a and M b ⊆ " M b " , this same as before Pareto front containing only the same mixtures of feasible molecules (e.g. ( M ka,2 , M kb,2,1 ) etc. or more generally of the form ( M ka,i , M kb,i , j ) ) would also be the only front satisfying condition (4). However, in this case mixtures of the form ( M ka,i ," M kb,i , j ") could never be identified hence completely missing feasible molecules of the form M ka,i , with the performanceassociated repercussions previously described. This indicates that the proposed approach serves to identify a set of high performance molecules in M a ,opt regardless of the structure of the molecules in M b ,opt and avoids the premature exclusion of molecules from M a ,opt that may lead to high performance mixtures. On the other hand, it might be argued that since structures " M kb,1,1 "," M kb,1,2 " are infeasible, mixtures with these molecules would be impractical. To address this argument it is required to consider potential ACS Paragon Plus Environment

15

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 61

individual Pareto fronts resulting from implementation of CAMD independently for each molecule included in M a ,opt . These Pareto fronts are conceptually illustrated in Figure 2b, indicating that for each molecule of the form M ka,i some of the corresponding molecules of the form " M kb,i , j " or M kb,i , j may be part of the overall Pareto front, while some others will only be optimum in the context of the individual Pareto fronts. For example, it appears that molecules of the form M kb,1,3 , M kb,1,4 are feasible molecular structures participating in mixtures with M ka,1 of nearly optimum performance, closely to the performance obtained by mixtures containing " M kb,1,1 "," M kb,1,2 " and belonging to the overall Pareto front, hence they can be used to replace the infeasible options within those mixtures. This is an argument in support of the necessity to capture molecules of the form M ka,1 . It should be noted here that the regions between the individual Pareto fronts and the overall Pareto front contain feasible solutions. There are points not shown in Figure 2b that dominate points of the individual Pareto fronts without belonging to the overall Pareto front and these points may lie in the gap between the individual and the overall Pareto front.

Stage 2: Mixture design This second design iteration serves to identify the optimum feasible molecular structures in M b ,opt and to also determine the optimum mixture concentration y1opt . This is necessary in case the Pareto front obtained in Stage 1 contains only infeasible structures in " M b ,opt " . In case that some of the structures in " M b ,opt " are feasible, Stage 2 can be avoided as the Pareto front already contains both optimum and feasible mixtures. However, it might still be necessary to implement Stage 2 for practical reasons. For example, feasible structures in " M b ,opt " might not be readily available for application (e.g. novel or commercially unavailable molecules), hence a second design iteration will serve to avoid them. The existence from Stage 1 of structural information resulting in high mixture performance can also be exploited to focus the search towards desired structures, avoid low performance areas and facilitate computations. Furthermore, as property predictions using GC methods involve uncertainty, a second 16 ACS Paragon Plus Environment

Page 17 of 61

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

design iteration may serve to review the molecules contained in M a ,opt and introduce higher accuracy in their expected mixture behavior by exploiting potentially available experimental data or other practical constraints. Stage 2 is approached by re-instating the chemical feasibility constraints in M kb (i.e. enforcing the satisfaction of conditions in equation (3)) and solving the mixture CAMD problem for each molecule of

{

}

the form M ka ,opt with optimization decision parameters d = mkb , Akb , y1 (Figure 3). There is no need to re-synthesize the first component of the mixture through vector M ka since molecules contained in M a ,opt have already been evaluated in terms of potential interactions with M kb in a much wider design space than in Stage 2 due to the lack of feasibility constraints. This means that no feasible structures in M kb evaluated with M ka during the optimization search in Stage 1 are excluded from this stage. Apparently, instead of repeating separately the CAMD design procedure for each molecule in M a ,opt , all molecules can be considered simultaneously as a finite set of discrete design options, hence implementing CAMD

{

}

only once with d = M ka ,opt , mkb , Akb , y1 . In both cases the design problem can be reduced by implementing upper and lower design limits based on the structural and property design insights derived using the information already available from Stage 1. In other words, regardless of the size of the obtained Pareto front in Stage 1, the optimization search in Stage 2 can be focused in promising areas of the design space. In this respect, Stage 1 is only concerned with the efficient and transparent identification of highly performing molecular structures in M a ,opt . The accurate identification of the optimum mixture concentration and the second component is left to stage 2, where the effort can be reduced as the user is allowed to review, interpret and analyze rich intermediate insights prior to exploiting meaningful conclusions.

Sensitivity analysis approach The mixtures obtained from the proposed design method will be required to satisfy practical process performance specifications. Uncertainty introduced in the design method due to approximations ACS Paragon Plus Environment

17

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 18 of 61

involved in the employed GC, thermodynamic or process models may result in over- or underestimation of the predicted thermodynamic or process behavior of the designed mixtures. Potential changes in the parameters employed for the calculation of mixture or process features (e.g. use of models of higher or lower accuracy at different design iterations, experimental data in place of GC property predictions etc.) are expected to result in larger or smaller mixture performance deviations. While several mixtures might be less sensitive to such changes, others are expected to significantly deviate from their calculated performance, eventually failing to meet the desired specifications. A systemic method is therefore required to classify the Pareto optimum mixtures in terms of their sensitivity to variations in thermodynamic property and process parameter estimates. This is addressed through a method facilitating the identification of parameters with high influence in the overall mixtureprocess system performance, the quantification of the overall system sensitivity with respect to these parameters and the incorporation of sensitivity metrics during the decision-making involved in the optimum mixture selection. The proposed sensitivity analysis method was originally developed44 for post-optimality analysis in nonlinear programming and subsequently used for the investigation of process systems controllability in view of external disturbances applied in reactive separation45 and separation systems46. The method was further applied to evaluate the controllability performance of extractive distillation processes utilizing different pure solvents47 derived from a CAMD-based method. In this work, the method is adapted to evaluate the sensitivity of mixtures in the presence of uncertainty in the employed thermodynamic property prediction and process models. It is implemented by imposing parametric variations on the prediction models around an optimal design point, d opt ∈ D opt , for the nominal parameter values. The aim is to calculate the impact of parameter variations on the employed objective functions. As Dopt represents the Pareto optimum mixtures, the following problem is solved for each mixture:

ACS Paragon Plus Environment

18

Page 19 of 61

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

Calculate F1 ( x, d , ε ) ,..., FNof ( x, d , ε ) h ( x, d , ε ) = 0

s.t.

(6)

g ( x, d , ε ) ≤ 0

{

x L ≤ x ≤ xU , d = d opt ∈ D opt ≡ M a ,opt , M b,opt , y1opt

}

where ε is the vector of parameters (GC, thermodynamic, process) on which variations are imposed with εnom corresponding to values of parameters at dopt and the entries of ε being scaled based on εnom to avoid biases. Note that the expression “calculate” in (6) refers to simulations for all realizations of vectors ε and Dopt and not to the implementation of optimization. This enables the calculation of the sensitivity matrix as follows:

 ∂ ln F1,1 ∂ ln F1, N p L  ∂ ln ε N p ,1  ∂ ln ε1,1  P= M O M  ∂ ln F ∂ ln FNof , N p N of ,1  L  ∂ ln ε1, Nof ∂ ln ε N p , Nof 

        

(7)

where d ln Fj ,i = dFj ,i / Fj ,i (j=1, …Nof ; i=1,…Np) represents a scaled transformation of the objective function values, while d ln ε i , j = d ε i , j / ε i , j is the scaled transformation of ε i , j . Sensitivity matrix P incorporates the scaled derivatives of the objective function values with respect to model parameters and constitutes a measure of the variation of the process model under the influence of infinitesimal changes imposed

on

(

model

parameters.

In

such

a

case

it

can

be

assumed

that

)

nom dF j ,i ≈ F j ,i x, d , ε inom , j + d ε i , j − F j ,i ( x, d , ε i , j ) with sufficiently small d ε i , j . Local sensitivity information

as depicted by matrix P might be useful for characterizing the variability of the obtained design results. However, such an approach will fail to capture nonlinear effects due to large in magnitude parameter variations and the interactions due to multiple simultaneous parameter variations. Furthermore, the existence of process constraints may lead to significant behavioral pattern changes for large multiple parameter variations. Therefore, an approach that allows for the investigation of the effect large and

ACS Paragon Plus Environment

19

Industrial & Engineering Chemistry Research

multiple parameter variations have on the optimal design and the corresponding objective functions is adopted. The local sensitivity matrix P is decomposed into major directions of variability by calculating the

Θ

eigenstructure of matrix PTP. The resulting eigenvectors

{ i }i=1 Np

are ranked ordered based on the

Θ magnitude of the corresponding eigenvalue. Eigenvector

1

that is associated with the largest in

magnitude eigenvalue represents the dominant direction of variability for the system, causing the largest change in the objective functions in a least squares sense45. A similar conclusion is derived for the second largest eigenvalue and so forth. The entries in the dominant eigenvector determine the major direction of variability in the multiparametric space and indicate the impact of each parameter in this direction. Having identified this direction it is not necessary to explore all directions of variability arbitrarily (i.e. combinations of parameters) hence reducing the dimensionality of the sensitivity analysis problem. Additional eigenvectors can also be considered in the sensitivity analysis capturing different directions of variability if a significant amount of the overall P matrix variability is described by them. Under the reasonable assumption that the eigenstructure of the local sensitivity matrix does not change considerably with the variation of parameters45, the following sensitivity metric is obtained for finite parameter variations: Np  Nof  F j ( x, d , ε (ζ ) ) − F j ( x, d , ε nom )   Calculate Ω(ζ ) = ∑ wiΩ  ∑    F j ( x, d , ε nom ) i =1  j =1   

s.t.

h ( x, d , ε (ζ ) ) = 0 g ( x, d , ε (ζ ) ) ≤ 0

ε (ζ ) − ε ε nom

nom



(8)

Θi

ζ =0

{

x L ≤ x ≤ xU , d = d opt ∈ D opt ≡ M a ,opt , M b ,opt , y1opt

}

Ω(ζ ) is called the sensitivity index, the fraction within the sum represents the relative change of objective function Fj , and ζ is a variation magnitude parameter indicating the range of the imposed

Θ

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 20 of 61

change in the direction of

. Weighting factor wiΩ introduces the significance of each objective

ACS Paragon Plus Environment

20

Page 21 of 61

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

function in the parametric sensitivity. The maximum change along the direction determined by each eigenvalue is dictated by the final value ζ lim of the ζ coordinate with ζ ∈ [ −ζ lim , ζ lim ] . This is usually selected based on process specific conditions regarding the expected parameter variability. The orthogonality condition for the system eigenvectors indicates that each eigenvector direction corresponds to a unique system variation mode. Investigation of the system behavior along more than one eigenvector directions enables the exploration of all significant underlying modes of variability for the system under study.

Implementation ORC modeling An ORC process consists of an evaporator (Ev) using a pure or mixed working fluid (wf) to draw heat from a heat carrier (hc). The fluid evaporates at high pressure and expands in a turbine to produce work. The low pressure vapor is then condensed (Cd) and pumped (Pm) back to the evaporator at high pressure to complete the cycle. The previously detailed16 ORC model is used in the course of CAMD to evaluate the system performance for different working fluid mixtures in view of the heat transfer limitations imposed by the their interactions with the heat carrier. Non-azeotropic mixtures exhibit a non-isothermal temperature profile during both evaporation and condensation, resulting in different ORC ORC values for the bubble and dew pressure. For the estimation of the maximum Pmax and minimum Pmin

cycle pressures, it is assumed that the mixture is in the state of saturated vapor at the evaporator outlet. ORC Ev , wf and the dew pressure Pdew |T Ev ,wf of the working fluid at the evaporator outlet As a result, Pmax out

temperature can be considered equal. In the condenser it is assumed that the working fluid is in the state ORC of saturated liquid at the outlet. Therefore Pmin is considered equal to the bubble point pressure

Cd , wf Pbub |T Cd ,wf of the working fluid at the temperature of the condenser outlet. The estimation of out

Ev , wf Cd , wf Pdew |T Ev ,wf and Pbub |T Cd ,wf is performed through vapor-liquid equilibrium (VLE) calculations using an out

out

Ev , wf equation of state. The temperature Tout of the evaporator outlet is calculated as follows:

ACS Paragon Plus Environment

21

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

Ev , wf Tout = TinEv ,hc − ∆Tmin

Page 22 of 61

(9)

where TinEv ,hc is the inlet temperature of the heat carrier and ∆Tmin a minimum desired temperature Cd , wf difference in the heat exchanger. The temperature Tout of the condenser outlet is set to a pre-specified

value, together with TinEv ,hc and the mass flow rate m hcf of the heat carrier. The mass flow rate of the working fluid is evaluated from a material and energy balance in the evaporator. In addition, it is assumed that a pinch point occurs in the evaporator when the working fluid starts to boil; therefore the Ev , wf temperature T Ev ,hc |T Ev ,wf of the heat carrier corresponding to the bubble point temperature Tbub of the bub

working fluid in the evaporator is calculated as follows: Ev , wf T Ev ,hc |T Ev ,wf = Tbub + ∆Tmin

(10)

bub

The turbine performance is evaluated under the assumption that the isentropic and mechanical efficiencies are constant and independent of the employed working fluid. The produced work W Tr is calculated as follows:

(

Tr wf W Tr = m wf ⋅ H out |( P ,T )Tr − H inwf |( P ,T )Ev f ⋅η out

out

)

(11)

wf where H out and H inwf are the enthalpies of the outlet and inlet streams calculated at the corresponding

temperature and pressure conditions noted in the equation and η Tr is the isentropic efficiency of the turbine. The work required by the pump is evaluated by the following simple expression: ORC ORC W Pm = m wff ⋅ (1/ ρ wf ) ⋅ ( Pmax − Pmin ) / η Pm

(12)

where ρwf is the average liquid density of the working fluid between the pump inlet and outlet and η Pm is the pump efficiency assumed equal for all working fluids. An enthalpy balance could also be used instead for the calculation of W Pm . The heat removed in the condenser and absorbed in the evaporator is calculated from the following equations:

(

wf Q Ev = m wf H out |( P ,T )Ev − H inwf |( P ,T )Ev f out

in

) ACS Paragon Plus Environment

(13)

22

Page 23 of 61

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

(

wf Q Cd = m wf H out |( P ,T )Cd − H inwf |( P ,T )Cd f out

in

)

(14)

Mixture design optimization problem formulation The multi-objective optimization CAMD approach employed in this work requires the utilization of a set of objective functions to be used during the screening and design stages. The consideration of an ORC model in the course of CAMD enables the utilization of process-related properties as objective functions together with mixture and pure-component properties as additional functions and constraints. There are numerous properties that should be considered for the design and selection of working fluid mixtures for ORC processes. Several of them are considered in this work to illustrate the implementation of the proposed methodology, which is independent of the selected properties as it can handle any desired performance measure or constraint. The mixture design problem is formulated as follows:

max nTh , nEx

(15)

min

RF1 , RF2

(16)

s.t.

ORC < Tcmix Tmax

(17)

ORC > max {Tm ,1 , Tm ,2 } Tmin

(18)

Tr Tout > Tdew |PTr

(19)

ORC < Pcmix Pmax

(20)

out

ORC min

P

>1

∞ sat 1 1 T2, fp

γ P | sat 1, fp

P

(21) >1 ∧

∞ sat 2 2 T1, fp

γ P | P2,satfp

CH2, >CH-, >CCF2, >CF-}. A maximum of 16 groups are allowed in each molecule. Boiling point temperature Tb, critical temperature Tc and pressure Pc of the pure components are the properties varied in the sensitivity analysis stage (vector ε) in a maximum allowable range of [{Tb,Tc}±30 K] and [Pc±10 bar] for each component. The overall design problem is solved twice throughout from CAMD to sensitivity analysis using both the Joback and Reid (JR)55 and Constantinou and Gani (CG)56 GC methods to calculate the three pure component properties ACS Paragon Plus Environment

26

Page 27 of 61

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

of interest. This serves to consider the variability in property predictions during CAMD and obtain an inclusive set of mixtures, prior to implementing the proposed sensitivity analysis approach in optimum mixture selection. These methods are chosen as they are widely utilized covering a broad range of properties, however more recent updates57 can also be used instead. While the ranges employed in the sensitivity analysis reported in this work are chosen randomly to illustrate the effects of uncertainty, a detailed discussion on uncertainty quantification in GC-based predictions is presented in literature39. The method of Lee-Kesler58 is used for the necessary ORC calculations with Plocker mixing rules59 and appropriate corrections for polar fluids60, based on a published implementation61. The required pure component critical properties and boiling point are calculated as noted previously, the ideal gas molar heat capacity is calculated based on JR55 and the binary interaction parameters are calculated based on published correlations employing both analytical62 and empirical methods63. Generally, the developments proposed in this work are independent of the employed methods and models, so others can also be used instead. Several mixtures were tested and compared well with the Peng-Robinson equation of state considered through the commercial process simulation package ProMax®64.

Results and discussion Structural and property mixture characteristics The implementation of stage 1 results in 19 different molecules in Ma,opt, identified as part of the Pareto front (Table 1). Among these molecules, 13 are identified as a result of implementing both JR and CG methods, 4 are identified only by the JR method and 2 only by the CG method. This is an indication that the variability observed in the property predictions performed by two major methods affects the type of resulting molecular structures. On the other hand, the molecules shared as a result of the two methods are considerably more. Although a wide range of potential molecular structures is covered by the two methods, the performed predictions result largely in coinciding molecules. In principle, the effects of the observed variability may be systematically addressed during molecular design using an uncertainty-based method35. However, the use of the two GC methods is sufficient

ACS Paragon Plus Environment

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

Page 28 of 61

within the scope of this work to both obtain an inclusive set of potentially optimum ORC mixtures and illustrate the proposed design approach. Several molecules shown in Table 1 represent common structures available in data repositories (e.g. webbook.nist.gov) such as Perfluoropropane (1), 1,1,1,3,3,3-Hexafluoropropane (2), 1,1,1-Trifluorobutane (8), Neopentane (CG1) etc. All molecules are alkanes, fluorinated alkanes or fluoromethoxyalkanes with flammabilities ranging from weakly or non-flammable fully fluorinated alkanes to increasingly flammable common pure alkanes. Furthermore, some of the molecules reported in Table 1 are also captured in “Mb,opt” as feasible structures of the stage 1 Pareto front. For example, 1,1,1,3,3,3Hexafluoropropane (2) captured in Ma,opt forms mixtures with 1-Fluoromethoxy-1,1,2,2,2-pentafluoroethane (JR1) and 1-Fluoromethoxy-2,2,2-trifluoro-ethane (JR2), which are captured in “Mb,opt”. At the same time molecule (2) is also captured in “Mb,opt” forming a mixture with Perfluoropropane (1) of Ma,opt. Fluoromethoxy-methane (10) is also captured in “Mb,opt” forming a mixture with Perfluoropopane (1), while 1-Fluoromethoxy-trifluoromethane (3) is also observed to form a mixture with 1,1,1,2,2Pentafluoropropane (5). On the other hand, all molecules of Table 1 with moderate to low flammability characteristics (less than the RF1 of molecule 11) are also subsequently identified in the Pareto front of stage 2 as part of Mb,opt and in various combinations with molecules of Ma,opt. This is partly illustrated through Table 2 (detailed later) showing the selected final mixtures after implementation of stage 2 and the proposed sensitivity analysis. For example, 1-Fluoromethoxy-propane appears as part of Mb,opt and has also been obtained as part of Ma,opt in stage 1. By a similar token, 1,1,1-Trifluoro-propane in Ma,opt has also been obtained in Mb,opt at various mixtures although it was not selected in the final results shown in Table 2 due to performance reasons subsequently elaborated in Figure 2. This is evidence that molecules obtained in Ma,opt of stage 1 facilitate the appearance of highly performing mixtures in a Pareto sense hence they are repeated either in “Mb,opt” or after implementation of stage 2 in Mb,opt. However, the consideration of a wider set of structural groups in CAMD is required to further verify these findings.

ACS Paragon Plus Environment

28

Page 29 of 61

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

It should be noted that the effects of variability due to the employed GC property methods are also verified from the molecular structures obtained in “Mb,opt”. Figure 4 shows the distribution of the groups contained in “Mb,opt” and obtained by both CG and JR methods. Structures obtained using the CG method indicate a higher frequency of “-CH3” and “>CH2” groups, while the JR method indicates a higher frequency of “>CH-”, “>C 1 ). This means that if the three investigated properties were

to receive values different than the ones corresponding to this ζ the mixture would result in an ORC system requiring vacuum operation in the condenser. For example, in the case of M6 for ζ=-0.08 the values of the three properties of 1,1,1-Trifluoropentane exactly before the pressure constraint is violated are [Tb2;Tc2;Pc2]= [330.83;480.24;27.36] with nominal values [Tb2;Tc2;Pc2]nom= [326.63;486.00;27.39]. Note that the imposed variations are still at ζ=-0.08 (i.e. close to the nominal points), however Tc2 and Tb2 are considerably different from Tcnom and Tbnom compared to the very small variation observed for Pc2 2 2 Tb 2 Tc 2 Pc 2 with respect to Pcnom ] = [-0.73;0.67;0.056] 2 . This is because the eigenvalue direction [ θ1 ;θ1 ; θ1

corresponding to Tc2 and Tb2 indicate that they have a very high impact on the ORC thermal and exergy efficiencies, while the ORC performance is only slightly affected by Pc2 (i.e. lowest absolute eigenvalue). Also note that only the first eigenvector Θ1 = [θ1ε i ]i =1p is considered in this case study, as it N

sufficiently describes the overall system variability. Apparently Tb2 and Tc2 are the two most influential properties due to the highest contribution in the eigenvalue direction. Due to the different signs Tb2 and Tc2 move with opposite trends. Table 3 reports the pure component properties ordered based on decreasing impact in ORC performance for each designed mixture. The considered variability space is both multi-dimensional and multi-directional due to multiple parameters evaluated simultaneously towards different directions. ACS Paragon Plus Environment

32

Page 33 of 61

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

Variations have not been imposed on 1,1,1,3,3,3-Hexafluoro-propane and Neopentane due to the use of experimental values instead of GC predictions. As a result, they were not considered in the calculation of the eigenvalue directions hence they are not reported in Table 3. In principle such variations could also be considered, as uncertainty is involved in experimental measurements. Although the proposed method allows for such an investigation, it is beyond the scope of this work. Except for the violation of process constraints, in some cases the profiles of Ω(ζ) are discontinued because the resulting solutions are practically infeasible. For example the profile of M6 for ζ=0.21 is discontinued because of a resulting negative accentric factor value. This may be avoided if a systemic investigation is performed for property combinations yielding feasible solutions. For ζ0. In all other cases there is violation of non-negative variables (e.g. accentric factor and mass flowrate) except for the case of M10 where VLE calculations are not converged within the specified solver tolerance. However, M10 already appears to be a very sensitive mixture, at least for ζ