Successive LP Approximation for Non-Convex Blending in MILP

Subscriber access provided by Kaohsiung Medical University

Process Systems Engineering

Successive LP Approximation for Non-Convex Blending in MILP Scheduling Optimization using Factors for Qualities in the Process Industry Jeffrey Dean Kelly, Brenno C Menezes, and Ignacio E. Grossmann Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.8b01093 • Publication Date (Web): 18 Jul 2018 Downloaded from http://pubs.acs.org on July 23, 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 56 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

Successive LP Approximation for Non-Convex Blending in MILP Scheduling Optimization using Factors for Qualities in the Process Industry Jeffrey D. Kelly,† Brenno C. Menezes, ‡, * Ignacio E. Grossmann‡ †

Industrial Algorithms Ltd., 15 St. Andrews Road, Toronto, Ontario, M1P 4C3, Canada.

‡

Department of Chemical Engineering, Carnegie Mellon University, Pittsburgh, United States.

ABSTRACT We develop a linear programming (LP) approach for nonlinear (NLP) blending of streams to approximate non-convex quality constraints by considering property variables as constants, parameters or coefficients of qualities that we call factors. In a blend-shop, these intensive properties of streams can be extended by multiplying the material flow carrying out these amounts of qualities. Our proposition augments equality balance constraints as essentially cuts of quality material flow for each property specification in a mixing point between feed sources and product sinks. In the LP factor formulation, the product blend quality is replaced by its property specification and variables of slacks and/or surpluses are included to close the balance; these are called factor-flows and well-known in industry as product giveaways. Examples highlight the usefulness of factors in successive substitution by correcting nonlinear blending deltas in mixed-

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 56

integer linear models (MILP) and to control product quality giveaways or premium specifications in blend-shops.

Keywords: blend scheduling, blend-shop, successive MILP, product giveaway, factor-flows.

1. Introduction Large-scale blend scheduling optimization problems found in the process industries are prohibitively expensive to solve in full space mixed-integer nonlinear programming (MINLP) approaches since the initial steps of the binary variable relaxation in a nonlinear problem (NLP) is highly degenerate (i.e., many non-differentiated solutions). To overcome this problem, the main strategy in the literature is to naturally decompose the quantity-logic-quality phenomena into a two-stage solution procedure.1-7 In these algorithms, the scheduling assignments in quantity-logic relationships are solved first in a mixed-integer linear problem (MILP) by neglecting the non-convex blending constraints commonly referred to as pooling.8 Then, a quantity-quality NLP problem is solved for fixed discrete or binary logic decisions found previously. The major drawback of this partition into MILP and NLP models is the possible gap between their objective functions. In addition, by setting up inconsistent discrete decisions, poor results and infeasibilities in the NLP solution may arise and is the reason we keep all logisticsfeasible solutions found during the MILP solution. Alternatives in the literature to improve the decomposed MILP sub-problems approximate the non-linearities from the NLP blending constraints in the mixed-integer solution with linear functions. This includes the following: a) a successive MILP engine by iteratively updating corrections on product quality bias (the difference between the linear and nonlinear calculated


2

Page 3 of 56 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


correlations of the blends defined before the optimization);9 b) piecewise McCormick envelopes to linearly under- and over-estimate bilinear terms;10-13 and c) multiparametric disaggregation to generate a tighter problem.2,14 However, the application of these approximations is limited due to the simplification of preferred blend recipes in a) and the considerable increase in the number of variables in b) and c). Although the tightness of the proposed factor model can be higher than the one found in the mentioned literature, the successive substitution procedure introduced later in this paper reduces iteratively the gap between the factor approach and the actual value of the blended quality stream. There are several other works covering the class of blending or pooling problems. Recent papers have shown that a linear approximation of this nonlinear programs is exact in some cases by using parametric structure of the objective function15 or that an MILP can generate an approximation with better performance than an NLP in other cases.16 Several works on the global optimization of non-convex quadratically-constrained optimization problems using piecewise-linear relaxations and other partitioning schemes can be found in the literature as well.17-21 Instead of neglecting quality information in the decomposed MILP sub-problems to avoid full space MINLP approaches in blend scheduling optimization problems, we propose a linear programming (LP) approximation for non-convex quality constraints from the blending of streams by using constants, parameters, coefficients, that we call factors. Factors can represent various types of qualities such as densities, components, properties and conditions including property blend indices or numbers. This approach is useful in blend-shops or where any controlled mixing in a process network takes place, whereby these intensive qualities of streams or factors are extended by multiplying them by the material flow and performing a factor-flow balance around the point or equipment of blending (or blender). For each property between raw


3


Page 4 of 56

material or feed sources and product sinks, an equality constraint or cut of quality material flow augments the linear formulas. To guarantee an LP balance of the amounts of qualities of each property, minimum and maximum specifications on the products substitute the quality of the blended material and variables for slacks and surpluses with unitary coefficients (for less than and greater than constraints, respectively) are included in the model. Furthermore, these factor balances also help to fine-tune or update the recipes from the planning optimization downloaded to the scheduling optimization as lower (hard), upper (hard) and target (soft) bounds. Since planned recipes are known to be gross estimates, factors can play a vital role to better manage the blend recipe uncertainties using more timely process quantity and quality data i.e., updated information within the planning time-period. When applying or proxying LP factors for quality in MILP problems, this can drastically reduce the MILP-NLP gap in the two-stage or phenomenological decomposition5 as more quality information is programmatically transferred to the logistics optimization. If this surrogate formulation of quality information (LP factors) is neglected in the MILP problem, many degenerate logistics solutions may occur with very close or similar MILP objective functions for different sets of quantity and logic decisions such as source to sink assignments and flow sizings. This unfortunately will increase the time to reduce the MILP-NLP differences or gaps when some iterative or sequential strategy is implemented, as found in our crude-oil blending and scheduling optimization problem by updating yields of distillates in the unit-operation out-portstates from the NLP into the new MILP iteration.7 A similar proposition for calculating overall constituent flow balances has been applied to full topology water network systems to significantly improve the strength of the lower bound in global optimum models,22,23 although this is considered as an NLP by the bilinear terms of


4

Page 5 of 56 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


blending of water stream flows with different contaminant concentrations. In our work, the concentrations are constants, parameters or coefficients (factors) instead of variables. Therefore, this approach is valid only for the mixing point between feeds and products, whereby properties or other intensive variables are considered as fixed values (factors), and therefore do not propagate throughout the network or flowsheet. Although analogous LP reformulations are common in NLP production planning problems by simplifying or aggregating the real topology of process networks to a sole point of mixing with constant qualities, it can be included in scheduling problems under certain circumstances as in the feed and product blend-shops with pseudo-constant concentration of the blending components as we propose. The LP factor approach introduced in this work can be applied to any process industry blending problem such as in the preparation of material ingredients in the food and beverage industry, feeding of diverse metal concentrates prior to smelting in the metals processing industry, and in the production of fuel, lube, asphalt and petrochemical streams in petroleum refineries. The examples highlight the usefulness of factors in successive substitution by correcting nonlinear blending deltas in blend scheduling MILPs, to decrease the MILP-NLP gap, to control and optimize product quality giveaways, and to reduce binary or discrete degeneracy in the logistics models.

2. LP approximation to non-convex NLP blending All blend-shops with varying quantities (extensive amounts) and qualities (intensive properties, concentrations, etc.) of raw component blending materials, gives rise to a non-convex NLP

problem in the calculation of volume- and/or mass-based properties (i.e., and , respectively) due to the plus and minus multi-linear terms in the quantity and quality conservation balances.


5


Page 6 of 56

As an alternative to these NLP blending equations, we propose an LP approach using factors for the quality variables to be applied in a decomposed MILP + NLP blend scheduling solution strategy. The topic of blending and scheduling of multiperiod problems is widely covered in the recent process industry literature,9,24-30 although we first introduce the use of the factor-flows or giveaways of qualities in such problems. The novelty includes an effective successive MILP approach using the optimized factor-flows to calculate linear-nonlinear property deltas in a postoptimization approach instead of in a pre-optimization stage using trial-and-error initial recipes to converge as in Mendez et al.9 Before the addressing of the LP factor modeling, we first highlight the traditional LP reformulation using the illustrative example of a simple blend-shop as shown in Figure 1 considering , and , as volume flows of inlet raw materials and outlet products ,

respectively. In this small heavy fuel oil (HFO) blending flowsheet, the diamond structures, shapes or objects (◇) are the sources and sinks known as perimeter unit-operations, whereby LAGO, LCGO and VR represent physical raw or component material flows of feedstock standing for light atmospheric gasoil, light cycle gasoil, and vacuum residue, respectively. LAGO and LCGO are considered cutting-stocks or fluxes of light diluent to primarily reduce the viscosity of the VR for the HFO production. The rectangular shape with the cross-hairs (⊠�) representing the HFO blender (HFOB) in Figure 1 is a continuous-process unit-operation with a blackbox subtype. The arrow external-streams or paths (→), also known as routes, connections, lineups, transfers or movements.


6

Page 7 of 56 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


Figure 1. Heavy fuel oil (HFO) blend-shop.

2.1. Traditional LP reformulation In Figure 1, the general in-out (or - ) material balance in the blender HFOB is calculated

as ∑ , , . The main qualities of the blended materials or products leaving the HFOB

output are determined by the volume-based constraint 1 and the mass-based constraint 2, where SG refers to the specific gravity or density property. These complete or proper quality balances of volume- and mass-based variables in the blender’s outlet (,, and ,, , respectively)

represent bilinear and trilinear terms given the product of extensive amounts in volume (, ) and intensive properties (,, and ,, ) of stock or material flows.31 The other nonlinear term in eq

1 below is related to the product of the blended material’s specific gravity property ,, and the sum of the inlet arrow flows , , represented by , . In eq 2, the other nonlinearity (trilinear) is in the product of mass-based property ,, , and the sum of the inlet accumulated mass

,, , (product of SG and volume flow), yielding ,, , .

The common linearization practice, shown in the inequality constraints 3 and 4, considers the feed or component properties as fixed or constant (̅,, and ,, ) and the properties of the

blended product are replaced by lower (̅,, and ,, ) and upper (̅,, and ,, ) bounds of

property specifications.3,5,6,9 In these linear constraints 3 and 4, the property variables of the


7


Page 8 of 56

blend cannot be exactly calculated since they are approximate or surrogate quality balances modeled as inequalities. This means there are slacks and / or surpluses in the quality calculation modeled as factor-flows (factor times flow), except when the constraint is binding or active; then the corresponding slack and / or surplus is zero. For less-than-or-equal (≤) constraints, which are nonbinding, some slack or unutilized amount of resource applies. For non-active greater-than-orequal (≥) constraints, some surplus or extra amount is being produced or utilized. These slacks or surpluses transform inequalities into equalities, which occurs internally in all LP solvers and is known as the standard LP form. ,, , ,, , ∀ , ∈ , 1

,, ,, , ,, ,, , ∀ , ∈ , 2

̅,, , ≤ ̅,, , ≤ ̅,, , ∀ , ∈ , 3

,, ̅,, , ≤ ,, ̅,, , ≤ ,, ̅,, , ∀ , ∈ , 4

We are proposing to explicitly have these over- or under-amounts of qualities (factor-flows of slacks or surpluses) modeled in an MILP optimization problem to improve the solution of the blending and scheduling decision-making automation.

2.2. Proposed LP approximation using factors for qualities The proposed LP approximation configures quality variables as fixed, invariant or constant parameters or coefficients of properties or factors which do not propagate over the flowsheet,


8

Page 9 of 56 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


network or superstructure. In Figure 2, HFO-DENS, HFO-SULF, HFO-FLSH and HFO-VISC perimeter unit-operations represent the factor-flow slacks (difference/distance to upper bound) and surpluses (difference/distance from lower bound) corresponding to the quality specifications of HFO.

Figure 2. Heavy fuel oil (HFO) blend-shop using factors. The proposed LP constraint 5 equalizes the extended quality amounts of property (factors multiplied by flows) around the blender unit-operation HFOB with respect to its factor-flows leaving the slack or surplus outlet % . The factor or coefficient in % is considered as unitary. Therefore, the value of the slack or surplus factor-flows & ,, represents, respectively, the insufficient or excess amounts of qualities for the LP factor flow of each property .

',, , ',, , ( & ,, ∀ % , , 5

In the left side of eq 5, for each property at time period , input stream amounts , in volume

with constant factors for qualities ',, counterbalance the total amount , of the blended ACS Paragon Plus Environment

9


Page 10 of 56

material factor or fixed property specification ',, minus slacks or plus surpluses of the factorflow variables & ,, , in the right side of eq 5. In such a way, for an upper bound of property

specification of the blended material represented by its factor ',, , there can be a slack or

negative factor-flow variable & ,, to complete the balance in eq 5, and hence & ,, ≤ 0. For a

lower bounded property of the blender’s factor ',, , a positive factor-flow or surplus is needed

for such factor balance (& ,, ≥ 0). As we will observe in example 1 (Section 5.1), when these factors are calculated as indexed of properties for viscosity, flash point and Reid vapor pressure (RVP), the change of scale from property to property index may transform positive numbers to

negative. Therefore, in the case of property indices, the factor-flow is modeled as & ,, ≤ 0 and

& ,, ≥ 0 to avoid infeasibilities.

Factors are similar to other quality phenomena variables such as densities, components, conditions and other physical-property coefficients except that, as mentioned, factors are constant. Factors can also be used to model simple enthalpy and entropy balances found in steam utility network optimizations where the factors are interfaced or linked to rigorous process simulators to update the enthalpy and entropy as temperatures and pressures change. These invariant coefficients ̅,, or ',, of upstream unit-operation out-port-states can be different giving values in each time window of the discrete-time step formulation proposed in this paper. In the HFO blend-shop example, there are four properties: specific gravity, sulfur concentration,

flash point and viscosity. For , in volume flow, the factor of specific gravity is its own value since it is a volume-based property. The sulfur concentration factor is the product of specific gravity and sulfur concentration (since it is a mass-based property), and the factor of flash point and viscosity are represented by their property index (also considered volume-based). The


10

Page 11 of 56 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


correlations used to transform property to its index are detailed in the HFO problem configuration in the Supporting Information (lines 254-280).

2.3. LP factor versus NLP blending approximations Considering the transformations from properties or indices to their respective factors, the LP proposition of this paper addresses the factor-flow balance (factor multiplied by flow) across the blender unit-operation by considering fixed properties or factors of intensive quality values such as specific gravity (SG) and sulfur concentration (SUL). Although eq 5 considers volume flows, a similar equation may be applied to any material basis such as mass flow. Commonly, volume flow is used in the models by the virtue of the widespread volumetric measurement of liquid amounts (flows and inventories) found in the process industry. In the mixture of raw material streams, for those properties varying linearly with the base flow, the blended property intensive value is one-dimensional and proportional to the extensions of the intensive value of the property streams. By using the LP factor model in linear mixtures, the factor-flows of slacks or surpluses (or giveaways) divided by the accumulated base flow matches precisely the distance between the lower and upper bounds of the quality specification and the calculated property (or index) at the amounts of the inlet material flows found in the LP factor approach (these can be checked in Figure 3 for , , , and & ,, in volume flow and in Figure 4 for ,, , ,, and ,& ,, in mass flow).

The mixing, pooling or blending relationship of volume-based properties (e.g., specific gravity SG) varies linearly with the volume flows of the inlet materials , (Figure 3, left side), although

for mass-based properties (e.g., sulfur concentration SUL) the factor SUL varies linearly with the


11


Page 12 of 56

mass flow ,, (Figure 4, right side). The equivalent equations to transform eq 1 to 5 (in volume flow) to a mass flow balance are represented in eq 6 to 10, considering -,, and -,, as factor

in mass flow for the inlets and outlet of the blender unit-operation. ,, ,, ,, ,, ∀ , ∈ , 6 ,,

,,

,,

,,

,,

,,

∀ , ∈ , 7

,, ,, ≤ ,, ,, ≤ ,, ,, ∀ , ∈ , 8 ̅,,

,, ̅,,

≤ ̅,,

,,

̅,,

≤ ̅,,

,, ∀ , ̅,,

∈ , 9

-,, ,, -,, ,, ( , & ,, ∀ % , , 10

To illustrate volume and mass flow variations with respect to their respective linear varying properties, Figures 3 and 4 use the HFO blend-shop example (in Figure 2) whereby it is considered fixed amounts of LCGO = 500.000 m3/d and VR = 800.000 m3/d and only one degree-of-freedom as the LAGO (light atmospheric gasoil) amount in volume and mass flows. These fixed amounts of LCGO (light cycle gasoil) and VR (vacuum residue) are the results of the blend-shop example for the HFO production found in Section 5.1 (Table 2) with LP factor result for LAGO in 237.619 m3/d, although LAGO is varied in Figures 3 and 4 from 0 to 1200 m3/d considering the fixed inlet properties of the LAGO, LCGO and VR components (̅ ,, , ̅,, or ,, ) and the specification of the HFO product in Table 1 (see Section 5.1).


12

Page 13 of 56 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


Figure 3. Specific gravity (SG) of the HFO blended product.

Figure 4. Sulfur concentration (SUL) of the HFO blended product. Figure 3 plots SG in both volume and mass flows considering a) SG in volume base (

) using

the true, real or actual nonlinear blending formula for volume-based properties in eq 1; b) SG calculated using the LP factor approach (─) in both eq 5 for volume flow and eq 10 for mass flow. For the flow in volume, in the left side of Figure 3, the factor ',, or ',, for SG is its own

value SG, reducing eq 5 to an in-out total mass flow conservation formula by the summation of the product of ',, and , for the inlets, and ',, and , plus or minus the factor-flow & ,, for the outlet of eq 5. For the flow in mass, in the right side of Figure 3, the factor -,,


13


Page 14 of 56

or -,, for SG is 1/SG, reducing eq 10 to an in-out total volume flow balance formula by the

summation of the product of -,2/, and ,, for the inlets, and -,2/, and ,, plus or minus the factor-flow ,& ,2/, for the outlet of eq 10.

Figure 4 plots SUL in both volume and mass flows considering a) SUL in mass base (

) using

the true nonlinear blending formula for mass-based properties in eq 2; b) SUL calculated using the LP factor approach (

) using both eq 5 for volume flow and eq 10 for mass flow. For the

flow in volume, in the left side of Figure 4, the factor ',, or ',, for SUL is the product given

by SG and SUL, reducing eq 5 to an in-out total mass flow conservation formula of SUL by the summation of the product of ',∙, and , for the inlets, and ',∙, and , plus or

minus the factor-flow & ,∙, for the outlet of eq 5. For the flow in mass, in the right side

of Figure 4, the factor -,, or -,, for SUL its own value SUL, reducing eq 10 to an in-out

mass flow conservation formula of SUL by the summation of the product of -,, and ,,

for the inlets, and -,, and ,, plus or minus the factor-flow ,& ,, for the outlet of eq 10. Comparing SG and SUL in volume and mass flows represented in Figures 3 and 4, respectively, the factor-flow approach accurately predicts the final or blended property quality for SG in volume flow (Figure 3, left side) and for SUL in mass flow (Figure 4, right side). However, for SG varying in mass flow (Figure 3, right side), there is an imperfection in the LP factor approach by using 1/SG as factor in eq 10 since the product of 1/SG and flow in mass gives a total volume flow, which is not truly conserved, as the summation of the inlet volumes of the blender maybe different from the outlet stream volume when liquids with different specific gravity (SG) are blended.


14

Page 15 of 56 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


For SUL varying in volume flow (Figure 4, left side), the LP factor approach uses the product SG and SUL as factor in eq 5 since the product of SG, SUL and flow in volume gives a mass flow of SUL such as in eq 2. Nevertheless, despite the summation of the inlet mass amounts of SUL around the blender is identical to the mass amounts of SUL in the outlet stream, because of mass flow of SUL is truly conserved, both SG and SUL are varying with volume flow. This is the main simplification in the LP factor model for mass-based properties with volume flow variations since the nonlinear curvature’s profile due to derivatives of bilinear and trilinear terms is neglected in the blended property variations, although this can be circumvented by an improved successive substitution with the support of the LP factor modeling as we introduce in the following section.

3. Successive substitution (SS) for nonlinear varying properties using factor-flows In previous successive LP approximations that replace NLP formulas for properties varying nonlinearly with the amounts of the components (such as for SUL in volume flow), the strategies for blending of streams use iterative property bias corrections initialized or pre-calculated before the start of the optimization.9 It demands a given initial recipe to compute the difference between the nonlinear formulas and linear volumetric blending, as well as an approximated specific gravity for mass-based properties. Instead, we calculate biases, corrections or deltas for any blending basis in a post-optimization stage by considering the LP factor results, specifically: a) Flows of components and blended product for the true NLP calculation using eqs 1 and 2. b) Factor-flows for the LP factor approach using eq 5. To illustrate the need of the SS method to properly compute sulfur concentration in a volume base flow modeling, we plot in Figure 5 the sulfur concentration of HFO product versus the


15


Page 16 of 56

amounts of LAGO considering the blend-shop example in Figure 2, where only LAGO volume in m3/d is considered as a degree-of-freedom. The complete data and results are found in the first example (Section 5.1) from which we extract insight to explain the usefulness of the LP factor theory to improve the successive LP approach that approximates the NLP blending correlations. The example in Figure 5 shows three iterations of the SS algorithm until convergence is achieved in both nonlinear calculation (to find on-specification properties) and in the prediction of property deltas (to avoid product giveaways after on-specification properties are reached). Although the method is explained using a single-period problem, the industrial-sized example demonstrates the scalability to a multiperiod case with hundreds of time-periods in the second example. The graphics in Figure 5 show the mass-based sulfur concentration property using the LP factor results to simulate the true nonlinear blended values of the products (

) by considering the

optimized flows in the LP factor approach and given properties of the raw materials in Table 1 (see Section 5.1). To successively correct the linear blended values of the products defined in the LP factor method (

), the SS iteration adds iteratively the property deltas, distances or shifts

at the virtual product specification until convergence is achieved. The property deltas are defined as the difference between the nonlinear calculated property (eqs 1 and 2) and the LP factor property (i.e., specification plus the factor-flows of slacks or surpluses divided by the accumulated basis, see the left side of eq 5). For properties and property indices governed linearly by volume-based relationships, the accumulated basis is the volume of the blended product. If the blended property is based in mass formulas such as for sulfur concentration and acidity, the accumulated basis is the volume multiplied by the specific gravity of the blended product (accumulated mass).


16

Page 17 of 56 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


Figure 5. SS method applied in the HFO blend for sulfur concentration (see Section 5.1). In the example, for the sake of the method illustration, the specification of sulfur concentration in the HFO product is arbitrarily reduced from its initial specification (a) at 3.25% to specification (b) at 3.12%, both % in mass. Considering the LP approach, in both specifications (a) and (b), the optimized flows of LAGO, LCGO and VR raw materials are the same, although the factorflows of sulfur concentration have reduced since the new specification (b) is closer to the LP factor results at the initial solution, equal to 3.0682% considering the reverse calculation of the factors from eq 5. Therefore, with the new specification (b), the blend is considered feasible at 237.619 m3/d of LAGO, the same results found at specification (a), indicated in Figure 5 as the initial solution. However, even if the LP factor solution is feasible in both specification (a) and (b) in the proposed modeling, there are under- and over-estimations of nonlinear properties of the blended material that can lead to infeasibilities or giveaways of quality. These inconsistencies are remarked in the HFO illustrative example for specification (b) at 3.12% as the true or real nonlinear simulation or calculation of the optimized flows (237.619 m3/d of LAGO and 500.000


17


Page 18 of 56

m3/d of LCGO) gives sulfur concentration at 3.1475% for SUL in mass. In this circumstance, the true quality calculation considering the LP factor solution would be infeasible. In the proposed successive-substitution (SS) algorithm, the infeasible blended value from the true simulated check (after an LP factor iteration) is reduced to substitute the uncertain LP factor blended property by new corrected values considering the property deltas from the past SS iteration. The example also explores the case of, after reaching a feasible solution (onspecification) in the nonlinear calculation with the use of the SS approach, the blended property yields a giveaway (over-specification) (see 2nd SS iteration in Figure 5, right side). In this case, the next SS iteration is performed by an interpolation between the current and last iteration. Figure 6 shows the SS algorithm considering the computation at each iteration of the delta corrections and the nonlinear calculation to check if the product is specified considering the results of raw material flows from the model using factors, both using an LP for blending or a MILP in blending and scheduling problems. In the automated SS search, new product specification with a correction of linear and nonlinear property delta is successively substituted in a new LP/MILP factor program by adding the current total delta from the past SS iteration until the product specification is reached as well as the additional delta of the last iteration converges lower than a tolerance (circa 0.0001).


18

Page 19 of 56 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


Figure 6. Successive-substitution MILP (SSMILP) algorithm using factors as quality proxies. By updating corrective deltas iteratively on the uncertain value of the blended material properties (or in the case the virtual specification given by spec’, spec’’ and spec’’’ in Figure 5), the improved successive substitution (SS) using factor-flows gives very close results to the nonlinear blending formulas as seen in the illustrative example. This procedure avoids infeasible solutions, checked at the nonlinear post-calculation using flows of raw materials and products, where the LP factor approach yields a feasible solution. Even if both the nonlinear and the linear postoptimization calculations give a feasible solution, the improved SS algorithm continues at each iteration to reduce the accumulation of the property deltas to near zero, therefore reducing the premiums or giveaways of qualities.


19


Page 20 of 56

As the successive substitution method starts the solution search at a point close to the accurate result by the addition of the factor quality cuts, the convergence of the method is guarantee considering that there are enough raw material quality properties to specify the product specifications. The nonlinear property delta between the true value (calculated in a postoptimization stage) and the one found by the LP factor approach is used, at each successive substitution iteration, to reduce the delta to lower than 0.0001. When the problem is only an LP, the speed of each iteration is expected to be fast, although for MILP problems it can be slower than without the use of the LP factor, since the quality factors actuate as cuts generating tighter problems than the original. Complex successive substitution to solve an NLP problem as a sequence of linear approaches can be found in Palacios-Gomez et al.32 and Kim and Ladson.33 The LP factor method is not guaranteed to global optimization. The addition of the factors or quality cuts approximates the NLP quality constraints as an LP approach, therefore it is useful in MILP problems that replaces MINLP methods, although an NLP solution can be performed after the MILP solution to recalculate the quantity and quality variables by fixing the resulting binary variable of the MILP. The successive substitution (SS) method based on the immutable factors for qualities of the inlet raw materials and the LP-NLP deltas, calculated in a post-optimization step, change the amounts of the raw materials by varying iteratively the virtual property (at the LP factor approach) of the resulted products. The method can fail if the initial point of the SS is far from the feasible space or if there are not enough raw material quality properties to specify the product specifications. However, as seen in the second example (Section 5.1), only with one SS iteration in the MILP problem the solution is near to the feasible space, reducing the possibilities of failure.


20

Page 21 of 56 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


4. Blend scheduling decomposition in MILP and NLP problems in UOPSS Figures 1 and 2 are reconstructed in the unit-operation-port-state superstructure (UOPSS)34,35 as shown in Figures 7 and 8, without and with the use of the LP factor, respectively. Besides, the unit-operations defined as perimeters (◇) for sources and sinks of materials, the unit-operations for continuous-process (⊠�) and the arrows or paths (→) connecting the units, there are circles with (⊗) and without () cross-hairs known as unit-operation outlet- and inlet-port-states, respectively. Port-state structures are unambiguous interfaces between up and downstream unit-

operations known as internal-streams, where movements can still be mixed in the in-port-states or split from the out-port-states . The index 6 represents out-port-states of upstream units

connected to the in-port-states , whereby 66 represents downstream in-port-states of other units connected to the sole out-port-state of the blender. In the UOPSS, the arcs can be attached to

the units in different time windows and points (or ports), and each port of attachment can still have a state (physical or not). In the very recent years, examples using UOPSS concepts found in Menezes et al.,36,37,43 Kelly et al.,38-40 Zyngier et al.41, Kelly and Menezes,42 and Franzoi et al.44 are solved for industrial cases, with many thousands of constraints and binary and continuous variables, within minutes using commercial solvers.

4.1. Blending formulas in UOPSS: NLP and LP reformulation using factors 4.1.1. NLP blending In Figure 7, the general in-out UOPSS material balance in the blender HFOB is calculated

as ∑ ∑7 7 ,, ∑ 77 , 77 , considering that multiple 7 ,, streams from upstream unit-operation

out-port-states 6 can be connected to diverse blender’s in-port-state . The output stream of a


21


Page 22 of 56

blender is normally dispatched for only one sink at a time, although it can be connected to more than one downstream in-port-states 66 .

Figure 7. Heavy fuel oil (HFO) blend-shop in UOPSS. The UOPSS volume-based constraint 11 and the mass-based constraint 12 in the blender’s out-

port-states (,, and ,, , respectively) represent bilinear and trilinear terms by the product of extensive amounts (7 ,, ) and intensive properties (7 ,, and 7 ,, ) of stock flows. These

streams leave out-port-state 6 of upstream unit-operations and arrive in the blender’s in-port-

states with ,, and ,, as property variables. If there is no mixing in these independent

inlets before they are blended in the HFOB unit-operation, the properties in the in-port-states are the same of the connected streams arriving from 6 ; therefore ,, 7 ,, and ,,

7 ,, . The general case for blending of streams considering in-port-states as mixers will be seen in the Section 4.2.3. The other nonlinear term in eq 11 is related to the product of the blended material specific gravity property ,, and each term of the sum of the arrow flows , 77 ,

leaving from the blender’s out-port-state . In eq 12, the trilinear terms are in the product of

mass-based property ,, and each term in the sum of their accumulated mass given by ,, , 77 , .

,, , 77 , ,, 7 ,, ∀ , ∈ , 11 77

7


22

Page 23 of 56 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


,, ,, , 77 , ,, ,, 7 ,, ∀ , ∈ , 12

77

7

4.1.2. LP factor blending The proposed LP reformulation in UOPSS in eq 13 equalizes extended quality amounts of

property (factors multiplied by flows) around the blender unit-operation HFOB with respect to its factor-flows leaving the slack or surplus out-port-states % (ohfo-dens, ohfo-sulf, ohfo-flsh,

ohfo-visc) in Figure 8. The factor in % is considered as unitary. Therefore, the value of the slack

or surplus factor-flows & ,, represents the insufficient or excess of amounts of qualities for the LP factor flow of each property respectively.

',, 7 ,, ',, , 77 , ( & ,, ∀ % , , 13

7

77

For each property to be reformulated in the UOPSS LP factor approach over a blender unit-

operation, raw material amounts ∑7 7 ,, arriving in different in-port-states with factors for

qualities ',, , in the left side of eq 13, counterbalance the total amount ∑ 77 , 77 , of the blended material factor or property specification ',, minus slacks or plus surpluses of the factor-flow

variables & ,, , the right side of eq 13. Although only one stream is connected to the blender

output in Figure 8, the summation ∑ 77 , 77 , in eq 13 represents several possible destinations of

the blended material as seen in the industrial example in Section 6 (see Figure 11).


23


Page 24 of 56

Figure 8. Heavy fuel oil (HFO) blend-shop using factors in UOPSS. For an upper bound of property specification in the product factor ',, , a slack or negative

factor-flow variable & ,, must close the balance in eq 13, so & ,, ≤ 0. A positive factor-flow

or surplus (& ,, ≥ 0) is necessary for a lower bounded property of the product blender’s factor ',, . To avoid infeasibilities, the factor-flow is modeled as & ,, ≤ 0 and & ,, ≥ 0 for

property indices as the transformation from property to property index may change positive numbers to negative.

4.2. Blend scheduling MILP and NLP problems With the statement of blend-shops using factors for qualities in the proposed LP approach around the mixing point of streams, we introduce the decomposition in MILP and NLP problems to replace, or more precisely, to approximate the solution of full space MINLPs for blend scheduling problems. The common part of the MILP and NLP problems considering binary or logic setups that are variable in the MILP and are fixed in the NLP is formulated in the following. The explanation of the blend scheduling optimization problems will be completed in


24

Page 25 of 56 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


the following sub-items: 4.2.1) MILP neglecting quality variables; 4.2.2) MILP using factor proxies or surrogates for qualities in LP blending formulas; and 4.2.3) NLP blending approach with fixed binary setups. Figure 9 shows a representative case of a blend-shop with sources of raw materials R1, R2 and R3, and their respective storage tanks S1, S2 and S3 connected to the same blender with modes of operation or grades A and B. Each mode produces grades of product with different properties or specifications. Tanks F1 to F4 hold the products to be delivered to the demand points D1 to D4. In addition to the objects in the HFO example of unit-operations, perimeters, port-states, and arrows (Figures 7 and 8), the blend scheduling problem must include storage vessels or tanks (∆) and their corresponding logic and logistics constraints.

Figure 9. Blend scheduling of fuels (blend-shops) in UOPSS. The objective function 14 below (suitable for both MILP and NLP problems) maximizes profit discounting the costs of feeds in R1 to R3 from the revenues of products (in the demand points D1 to D4), although in many cases the costs of the raw material streams are disregarded for


25


Page 26 of 56

producers of fuels since it is complicated to define reasonable costings of the intermediate streams dispatched downstream from the plant to the blend-shops. The UOPSS formulation, given by the objects and their connectivity as in Figure 9, are defined in eqs 15 to 30. In the

summations involving ports in eqs 18 to 29, 6 and 66 represent upstream and downstream portstates connected, respectively, to the in-port-state and out-port-state of unit-operations ,,

yielding a combination of 6 , , , 66 in-out flow of a unit-operation ,. The sets 89: , 8;< , 8= and 8> define unit-operations , for, respectively, raw materials, tanks, blenders and product

demands. The sets ? and @ define in- and out-port-states, respectively, and their indices are shown in the summations of eqs 18 to 29 for storage (A) and feed (B) tanks. The set @? represents the arrow flows connection an out-port-state to an in-port-state. For ∈ ℝD and E {0,1}:

8H I J K LMN, N, − LQRN, N, U 14 s.t.

N∈:O

N∈:ST

̅,, E,, ≤ ,, ≤ ̅,, E,, ∀ , ∈ @?, 15

̅N, EN, ≤ N, ≤ ̅N, EN, ∀ , ∉ 8;< , 16 Xℎ XXN, XXXN, EN, ≤ ℎN, ≤ ℎ EN, ∀ , ∈ 8;< , 17

7 ,, ≥ ̅ N, EN, ∀ , , ∈ 8> , 18

7 ∈Y&

7 ,, ≤ ̅ N, EN, ∀ , , ∈ 8> , 19

7 ∈Y&

, 77 , ≥ ̅N, EN, ∀ ,, ∈ 89: , 20

77 ∈Z[

, 77 , ≤ ̅N, EN, ∀ ,, ∈ 89: , 21

77 ∈Z[


26

Page 27 of 56 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


EN, ∀ ,, ∈ 8= , 22 , 77 , ≥ ̅ N,

77 ∈Z&

, 77 , ≤ ̅ N, EN, ∀ ,, ∈ 8= , 23

77 ∈Z&

7 ,, ≥ K̅, N, ∀ , , ∈ 8= , 24

7 ∈Y[

7 ,, ≤ K̅, N, ∀ , , ∈ 8= , 25

7 ∈Y[

, 77 , ≥ K̅, N, ∀ ,, ∈ 8= , 26

77 ∈Z&

, 77 , ≤ K̅, N, ∀ ,, ∈ 8= , 27

77 ∈Z&

ℎN, ℎN,\2 ( 7 ,, − , 77 , ∀ , ,, ∈ 8;< , 28 7 ∈Y]^

77 ∈Z]^

7 ,, , 77 , ∀ , ∈ 8= , 29

∈Z_` 7 ∈Y[

77 ∈Z&

N, , ,, , ℎN, ≥ 0; E,, , EN, {0,1} 30

The semi-continuous constraints 15 to 17 relate the binary and continuous variables of the flowsheet depicted in Figure 9 for: a) quantity-flows of the arrows ,, for , ∈ @?; b)

throughput-flows of the non-tank unit-operations N, (, ∉ 8;< ); c) tank-holdups or inventory

levels ℎN, (, ∈ 8;< ). If their binary variables E,, or EN, are active, these amounts (flows and holdups) vary between bounds: ̅,, and ̅,, for ,, , ̅ N, and ̅N, for N, (, ∉ 8;< )

XXXN, XXXN, and ℎ and ℎ for ℎN, (, ∈ 8;< ). In equations 18 to 23, when the binary variable EN,

of the unit-operation , (, ∉ 8;< ) is active, the summation of the quantity-flow of the arrows


27


Page 28 of 56

and arriving in or leaving from its port-states must be between the unit-operation bounds (̅ N,

̅N, ). Equations 18 and 19 represent the sum of the arrows arriving from the upstream out-port-

states 6 of the feed tanks ( 6 ∈ @% ) into the in-port-states of the demand points ( , ,) ∈ 8> .

Equations 20 and 21 represent the sum of the arrows leaving from the upstream out-port-states

of the raw material sources of ,, ∈ 89: into in-port-states 66 of the storage tanks ( 66 ∈ ? ).

Equations 22 and 23 are the sum of the arrows leaving from the out-port-states of the blender

(, ∈ 8= ) in direction to the downstream in-port-states 66 of the feed tanks ( 66 ∈ ?% ).

Blender unit-operations , (, ∈ 8= ) can have streams arriving in or leaving from independent

connected port-states. Equations 24 and 25 consider bounds of inverse yields (K̅, and K̅, ) and eqs 26 and 27 consider bounds of direct yields (K̅, and K̅, ) with respect to the unit-operation

throughputs N, of the blender. However, typically only one outlet flow is considered for

blenders, although they can be connected to several sinks, one operating at a time in a normal

blender operation. The quantity balance of inventory or holdup for unit-operations , of tanks (, ∈ 8;< ) is calculated in eq 28 considering the current holdup amount ℎN, of the heels or material left in the past time-period plus and minus the summation of the upstream and downstream connections to the tanks. For no accumulation of material in blenders, a material balance in eq 29 is necessary in these types of units. It should be mentioned that the UOPSS formulation only calculates material balances in unit-

operations N, (, ∉ 8;< ) with respect to the networked flows of the arrows ,, at bounds of inverse yields (K̅, and K̅, ) and yields (K̅, and K̅, ) as in eq 24 to 27. The setup operation

balances in eqs 18 to 23 are with respect to the possible flows in and out of a unit-operation , (, ∉ 8;< ) when its setup EN, is active. This can be considered a semi-continuous equation ACS Paragon Plus Environment

28

Page 29 of 56 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


linking the network of flows to their connected unit-operations instead the unit-operations to themselves as in eqs 16 and 17.

4.2.1. MILP blend scheduling without factors for qualities The MILP problem neglecting NLP quality constraints and without factors of qualities includes eqs 14 to 30 and the following logistics constraints in eqs 31 to 44 that involves: a) structural transition constraints and modes of operation selections; b) unit-operations’ temporal transitions, c) multi-use of ports, d) uptime or minimal time of unit-operations, and e) tanks in fill-draw delay operations (Zyngier and Kelly45). In eq 31, if the binary or setup variables of interconnected unit-operations ,6 and , are active

(,6 is the upstream unit-operation to ,), the binary variable E,, of the arrow connecting them within out-port-states 6 of ,6 and in-port-states of , is implicitly active. This aggregated

structural transition constraint is a collection of 4 objects (,6 , 6 , , ,) which is a logic valid cut that facilitates the tree search in branch-and-bound methods and can be disaggregated if required. In eq 32, for all physical units, at most one unit-operation , (as EN, for procedures, modes or

tasks) is permitted in bN at a time. The set bN represents the unit-operations , within the same

physical unit. EN7 , ( EN, ≥ 2E,, ∀ ,6 , , , ,, 31

EN, ≤ 1 ∀ 32

N∈c


29


Page 30 of 56

The temporal transition constraints 33 to 35 from Kelly and Zyngier46 control the operation of

the semi-continuous blender HFOB, where the setup or binary variable EN, (, ∈ 8= ) actives the start-up, switch-over-to-itself and shut-down variables (dReN, , dRN, and dRfN, ,

respectively) that are also binary variables. Equation 34 may also be disaggregated if required and eq 35 encourages the integrality of the temporal transition logic variables. EN, − EN,\2 − dReN, ( dRfN, 0 ∀ , ∈ 8= , 33 EN, ( EN,\2 − dReN, − dRfN, − 2dRN, 0 ∀ , ∈ 8= , 34 dReN, ( dRfN, ( dRN, ≤ 1 ∀ , ∈ 8= , 35

In the multi-use procedure in eqs 36 and 37, the lower and upper parameters bAg, and bAg,

coordinate the utilization of the out-port-states ( ∈ @h ) by their connected downstream in-portstates ′′. This occurs in the out-port-state of blenders to avoid the split of a mixture to different feed tanks at the same time. Equations 38 and 39 impose multiple utilization bounds bAg, and

bAg, in the in-port-states ( ∈ ?h ) with their connected upstream out-port-states 6 .

E, 77 , ≤ bAg, EN, ∀ ∈ @h , 36 77

E, 77 , ≥ bAg, EN, ∀ ∈ @h , 37 77

E7 ,, ≤ bAg, EN, ∀ ∈ ?h , 38 7

E7 ,, ≥ bAg, EN, ∀ ∈ ?h , 39 7


30

Page 31 of 56 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


Equations 40 to 42 model the run-length or uptime use considering bj and bj as the lower and upper bounds of time of using, respectively, and klm as the end of the time horizon with ∆t as time-step, where eq 42 is the unit-operation uptime temporal aggregation cut constraint and number of period is o .45,46

q; ` \2

dReN,\pp ≤ EN, ∀ , ∈ 8= , > 1 | − ≥ o 40

2

q; t ∆

EN, ≤

bj ∀ , ∈ 8= , < klm − bj 41 ∆

∆ dReN, ≤ o ∀ , ∈ 8= 42

The fill-draw delay constraints 43 and 44 controls, respectively, the minimum and maximum time between the last filling and the following drawing operations for the upstream 6 and

downstream 66 connections of a tank. The minimum and maximum fill-draw delays are ∆vNl and ∆vNwx , respectively, and the in- and out-port-states and are connected to the tank

(?@∈8;< ).39

E7 ,, ( E, 77 ,D ≤ 1 ∀ 6 , , , 66 such that ?@∈8;< , ∆>c

0. . ∆vNl , 1. . | ( < klm 43

E7 ,,\2 − E7 ,, − E, 77 ,\2 ≤ 0 ∀ 6 , , , 66 2

such that ?@∈8;< , 1. . − ∆vNwx | ( < klm 44


31


Page 32 of 56

4.2.2. MILP blend scheduling with factors for qualities The MILP problem using factors for qualities in LP blending formulas includes eq 13 for each property specification as the LP factor proposition of this paper, eqs 14 to 30 for the UOPSS formulation and eqs 31 to 44 for the logistics modeling. Here, instead of neglecting NLP blending constraints in the MILP problem, the LP reformulation using factors (eq 13) is applied as seen in Figure 10 for specific gravity (SG) and sulfur concentration (SUL).

Figure 10. Blend scheduling using factors.

4.2.3. NLP blend scheduling The NLP problem includes eqs 11 and 12 for volume- and mass-based properties in the out-portstate of blenders and the UOPSS formulation in eqs 14 to 30 for fixed binary setup variables.


32

Page 33 of 56 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


The NLP depends on the MILP solution to define the setup of unit-operations (including tanks) and unit-operation-port-state to unit-operation-port-state (the arrows), therefore only quantity and quality variables are determined.

The remaining NLP constraints for the UOPSS network include a) non-tank in-port-state as

mixers in eqs 45 and 46 and b) tanks , with holdups in eq 47 and outlet flows in eqs 48 and 49

(the equation for mass-based properties in tanks is skipped since it is similar to eq 47 by replacing by the product ).

,, 7 ,, 7 ,, 7 ,, ∀ ∉ ?;< , ∈ , 45 7

7

,, 7 ,, 7 ,, 7 ,, 7 ,, 7 ,, ∀ ∉ ?;< , ∈ , 46 7

7

N,, ℎN, N,,\2 ℎN,\2 ( 7 ,, 7 ,, − N,, , 77 , ∀ , ,, ∈ 8;< , 7

77

∈ , 47

,, N,, ∀ ,, ∈ 8;< , ∈ , 48 ,, N,, ∀ ,, ∈ 8;< , ∈ , 49

5. Illustrative Examples Before we start with the blend scheduling examples with and without the proposed LP approximation using factors for qualities, the illustrative case of a single-period blend-shop is presented to compare our LP proposition (Figure 8) with the NLP approach (Figure 7) for mixing point of streams. All examples are performed in an Intel Core i7 machine at 3.41 GHz (in 8 threads) with 64GB of RAM and use the structural-based UOPSS framework found in the


33


Page 34 of 56

semantic-oriented modeling and solving platform IMPL (Industrial Modeling and Programming Language) from Industrial Algorithms Limited.

5.1. Heavy fuel oil (HFO) blend-shop: LP using factors and NLP approaches Also known as Fuel Oil #6 and Marine Bunker C, Heavy Fuel Oil (HFO) is typically a blend of Light Atmospheric Gas Oil (LAGO) and Light Cycle Gas Oil (LCGO) with Vacuum Residue (VR) as seen in Figures 1 to 4. For oil-refineries without a Delayed Coker, Visbreaker, etc. there can be a significant amount of VR produced from a Vacuum Distillation Unit (VDU) which must be mixed with a flux or cutter-stock to reduce the viscosity of the residual fuel oil. HFO is considered to be a by-product since its market price is less than the cost of a crude-oil, and therefore it is important not to over-flux with LAGO/LCGO which can be disposed into more valuable refined products. For this HFO blending example in the proposed LP approximation with factors, we have configured four factors, DENS, SULF, FLSH and VISC as seen in Figure 8 (using the UOPSS framework). The out-port-states HFO-DENS, HFO-SULF, HFO-FLSH and HFO-VISC are hypothetical flow variables (factor-flows) representing the excess or insufficient property specifications and these flows have the quality bounds configured. Density blends linearly by volume and sulfur blends linearly by mass. However, flash-point and viscosity require what are known as blending indices which after the transformation blend linearly by volume (see the HFO flowsheet set and configuration in the Supporting Information). Table 1 shows the initial properties before they are transformed into factors. For this example, we have a flow of VR of 800.000 m3/d which must be blended in its entirety to meet required HFO quality specifications. The lower and upper flow bounds for LAGO are LCGO are 0.0 to


34

Page 35 of 56 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


500.0 m3/d respectively where their respective costs are $65.00 and $55.00 and the price for HFO is $40.00 (all per barrel). The measurement biases for all factors are normalized to zero (0.0) although these can be used to offset estimated process and measurement uncertainty in the system. Table 1. Example 1's data for heavy fuel oil (HFO) blending. DENS g/cm LAGO LCGO VR HFO

3

0.85 0.98 1.03 ≤ 1.01

SULF %w

FLSH ºC

VISC cSt at 50ºC

0.80 1.00 5.00 ≤ 3.25

50 78 200 ≥ 65

5 10 250000 ≤ 380

The values in bold in Table 2 are the resulting flows and factor-flows (FF) in the HFO production example using an LP with a single time-period. As expected, the solution maximizes the injection of LCGO ($55.00) given that it is cheaper than LAGO ($65.00) where LAGO is added to ultimately respect the upper viscosity bound of 380 cSt (mm2/s) at 50 degrees ºC (or 122 ºF). Given that the LP problem has two degrees-of-freedom i.e., LAGO and LCGO flows, the LP automatically determines that the optimal solution is to use the LCGO upper bound (500.000 m3/d) and the HFO-VISC as the active constraints or basis, using a minimal of 237.619 m3/d of LAGO to match the HFO viscosity. Therefore, the HFO production is equal to 1,537.619 (for VR fixed in 800.000 m3/d) with profit objective function of $116,735.80/d. Table 2. Factor-flows (FF) for the LP HFO blending (Flowsheet in Figure 8).


35


Flow m 3/d 237.619 500.000 800.000 1537.619

LAGO LCGO VR HFO Specification Slack or Surplus Property or Index Delta Calculated Property or Index (linear) Calculated Property or Index (true nonlinear) % (linear - nonlinear)/nonlinear

FF-DENS m 3/d . g/cm3 201.9761 490.0000 824.0000 1515.9761 1552.9952 -37.0190 -0.0241 0.9859 0.9859 0.00%

FF-SULF m 3/d . g/cm3 . %w 161.5809 490.0000 4120.0000 4771.5809 4845.3449 -73.7640 -0.0487 3.0713 3.1475 2.42%

Page 36 of 56

FF-FLSH m 3/d . [Index] 72678.4244 25333.4371 246.0289 98257.8903 172059.5622 -73801.6718 -47.9974 63.9026 63.9069 0.01%

FF-VISC m 3/d . [Index] -7456.9071 -14046.6800 -1777.5034 -23281.0904 -23281.0904 0.0000 0.0000 -15.1410 -15.1410 0.00%

The slacks for HFO-DENS, HFO-SULF and HFO-VISC found in the LP solution are, respectively, -37.0190, -275.6534 and 0.0. The value of the surplus for HFO-FLSH is 73,801.6718 considering this the index of the flash property. As the slack for HFO-VISC is zero, this is an indication that the viscosity is a limiting factor. Although the surplus of HFO-FLSH is a negative number due to the transformation from property to index of property, the flash point property surplus is positive. Considering the LP factor approach results in bold in Table 2 (flows and factor-flows) and the given properties in Table 1, the linear and nonlinear calculations of properties (or indices) are also compared in Table 2. As expected, the linear and nonlinear differences for those properties or indices varying linearly with volume is zero (as for density and indices of viscosity and flash), although for sulfur concentration there is an error around 2.42%. To further explore this example, we arbitrarily reduced the sulfur concentration from spec. (a) at 3.25% (see Table 1) to a new specification, spec. (b), at 3.12%. Even though, the varying flow of raw materials LAGO and LCGO would be the same at both spec. (a) and (b) since the LP factor approach gives a feasible solution at 3.0682%, albeit the factor-flows or slack of sulfur concentration increases from -275.6534 (at 3.25%) to -73.7640 (at 3.12%) as the latter gets


36

Page 37 of 56 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


closer to the specification. However, when the nonlinear simulation or calculation is checked at the new spec. (b), it will give an infeasible solution, since 3.1475% is bigger than the spec. (b) upper bound of 3.12%. For this case we apply the successive substitution (SS) algorithm in Section 3 (see the steps to find the final solution in Figure 5). The SS at the new specification 3.12% demands three iterations to both get a feasible solution (when simulating the optimized flows of raw materials of the LP results into the nonlinear blend) as well as to reduce the nonlinear and linear accumulated property deltas to lower than the tolerance of 0.0001 (see Figure 5, right side). After such SS algorithm, as shown in Table 3, the amount of LAGO increased from 237.619 to 258.780 m3/d to guarantee the minimal use of the lighter diluent (LAGO) at the maximum sulfur concentration at 3.12% in mass, at spec. (b), resulting in a total amount of 1558.780 m3/d and a profit of $113,639.98/d. In this case, considering the linear and nonlinear differences for sulfur concentration there is an error around 2.63%, however the true nonlinear property is specified at 3.12%. Table 3. Factor-flows (FF) for the LP HFO blending after the SS algorithm. Flow m 3/d 258.780 500.000 800.000 1558.780

LAGO LCGO VR HFO Specification Slack or Surplus Property or Index Delta Calculated Property or Index (linear) Calculated Property or Index (true nonlinear) % (linear - nonlinear)/nonlinear

FF-DENS m 3/d . g/cm3 219.9630 490.0000 824.0000 1533.9630 1574.3678 -40.4048 -0.0259 0.9841 0.9841 0.00%

FF-SULF m 3/d . g/cm3 . %w 175.9704 490.0000 4120.0000 4785.9704 4912.0275 -126.0571 -0.0822 3.0378 3.1200 2.63%

FF-FLSH m 3/d . [Index] 79150.7642 25333.4371 246.0289 104730.2301 174427.4820 -69697.2519 -44.7127 67.1873 67.1923 0.01%

FF-VISC m 3/d . [Index] -8120.9781 -14046.6800 -1777.5034 -23945.1615 -23601.4897 -343.6718 -0.2240 -15.3650 -15.3615 -0.02%

As a comparison, using eqs 45-49 in Section 4.2.3, the NLP approach of the HFO blend-shop in Figure 7 is solved as well. The amounts of LAGO and LCGO in this case are the same as found


37


Page 38 of 56

in the LP blending with factors at spec. (a) (see Table 2) and at spec. (b), as seen in Table 3. The successive substitution algorithm is applied in spec. (b) since the nonlinear simulated check is infeasible only in this case (see SUL in mass at the initial LAGO amount in Figure 5). Furthermore, as in the NLP model the properties or indices of the blended material (,, and

,, ) are obtained, it is possible to check if the factor-flows of slacks or surpluses (& ,, ) when divided by their accumulated base and added in the lower (̅,, and ,, ) or upper (̅,, and

,, ) bounds of property or index specifications, give the same optimized ,, and ,, as in

the NLP problem. After the check, in the case of spec. (b) for sulfur concentration, the one where the SS iteration in applied, the calculated density and indices of flash and viscosity using the flow and factor-flows in the LP problem with SS and the NLP results are identical, as seen in Table 4. Table 4. LAGO amounts and HFO properties (DENS, SULF) and indices (FLSH, VISC) in the SS LP factor and NLP models.

LAGO (SS LP Factor) LAGO (NLP) % (LP - NLP)/NLP

Flow m 3/d 258.780 258.783 0.0012%

DENS 3

kg/m 984.0792 984.0790 0.0000%

SULF %w 3.1200 3.1200 0.0000%

FLSH [Index] 67.1873 67.1923 0.0074%

VISC [Index] -15.3650 -15.3615 -0.0232%

5.2. Blend Scheduling in MILP with factors At this point we explore the blend scheduling problem in Figure 9. The phenomenological decomposition of the MINLP model uses the heuristic of solving pure MILP scheduling problems first, then to fix their binary results in the NLP quality sub-problems. However, rather than neglecting qualities in the MILP models, as in the decomposed blend scheduling examples


38

Page 39 of 56 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


without factors (BS) as in Figure 9, equality balances are included by using factor-flows in the in-port-states and out-port-states of the blender with factors (BSF) as seen in Figure 10. In the logistics of the blending and scheduling problems (BS and BSF), there is a minimum and maximum uptime or block-time for using of the blenders, meaning that once its use has started, it should last between 1 and 3 hours. Also, a fill draw delay of 1 hour is defined for F1 to F4 as the minimal time to permit the mixing or homogenization of the raw material after the last in-port movement from the blender HFOB. The holdups or inventory levels of all tanks (S1, S2, S3, F1 and F2) have lower and upper bounds of 10 and 200 m3. The quality problem has constant specific gravity (SG) and sulfur mass concentration (SUL) properties in the raw material sources (R1 to R3) and their respective storage tanks (S1 to S3) with the following values: R1 and S1 (SG = 0.800; S = 0.950; cost = $30); R2 and S2 (SG = 0.825; SUL = 1.000; cost = $20); R3 and S3 (SG = 0.880; SUL = 1.500; cost = $10). The same blender operates in modes A and B to produce, respectively, D1 (SG = 0.815; SUL = 0.970; price = $70) and D2 (SG = 0.885; SUL = 1.100; price = $60) products. Table 5 shows the results of the blend scheduling with factors (BSF) for 8 hours with 1h as timestep considering the sulfur concentration (SUL) in the MILP factor approach convergence after 5 successive substitution (SS) iterations. As a comparison, the MILP-NLP gap in the blend scheduling without factors (BS) is 24.49% for a MILP solution of $460.00 and NLP of $347.34 (not presented in Table 5), although the MILP-NLP gap is successively reduced in five successive substitution (SS) iterations to 0.03% when using factors instead of neglecting the qualities as in the BS problem. The MILP solvers use best commercial solvers CPLEX (12.71) and GUROBI (7.5.1) and the NLP solver is connected to them using the SLPQPE (Sequential


39


Page 40 of 56

Linear Programing Quadratic Programming Engine) technique found in IMPL. In all MILP cases in this example, the relaxation gap is set to 1%.

Table 5. Sulfur concentration (SUL) in scheduling with factors (BSF) at each iteration. Objective ($) st

1 SS nd 2 SS rd 3 SS th 4 SS th 5 SS

MILP 354.20 349.25 348.00 347.59 347.44

MILP-NLP Gap 1.98% 0.55% 0.19% 0.07% 0.03%

A SUL 0.9739 0.9714 0.9705 0.9702 0.9701

B SUL 1.1569 1.1044 1.1003 1.1000 1.1000

A

B SS

delta 0.0039 0.0014 0.0005 0.0002 0.0001

SS

delta 0.0569 0.0044 0.0003 0.0000 0.0000

A

B

delta 0.0039 0.0053 0.0059 0.0061 0.0061

delta 0.0569 0.0613 0.0616 0.0616 0.0616

Table 6 shows the SS iteration results of the inlet raw material flows from the storage tanks S1 to S3 for 8 hours with 1h as time-step. The product demands are considered open in all time windows and the preferred stream (6 time-windows) is the product from the blender in the mode of operation A (SG = 0.815; SUL = 0.970). This preferred product occurs since the 1st SS iteration without any correction in the virtual or uncertain property specification of sulfur concentration, although the amounts of each inlet stream vary in each blend A and B to specify the SUL in both products with the use of the SS algorithm.

Table 6. Inlet flows in the blender unit-operations (grades A and B).


40

Page 41 of 56 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


st

1 SS S1 S2 S3 nd

2 SS S1 S2 S3 rd

3 SS S1 S2 S3 th

4 SS S1 S2 S3 th

5 SS S1 S2 S3

1 A 53.00 47.00

2 B

A 57.87 42.13

70.00 30.00 A 57.87 42.13

A 59.68 40.32

A 59.68 40.32

A 60.34 39.66

A 60.34 39.66

A 60.59 39.41

A 60.59 39.41

3 A 53.00 47.00

4 A 53.00 47.00

B

A 57.87 42.13

80.17 19.83 A 59.68 40.32

B 81.01 18.99 B 81.01 18.99

B 80.95 19.05 A 60.34 39.66

A 60.59 39.41

5 B 70.00 30.00 A 57.87 42.13

A 59.68 40.32 A 60.34 39.66 B 81.01 18.99

6 A 53.00 47.00

7 A 53.00 47.00

8 A 53.00 47.00

A 57.87 42.13

B

A 57.87 42.13

B 80.95 19.05 B 81.01 18.99 A 60.59 39.41

80.17 19.83 A 59.68 40.32

A 59.68 40.32

A 60.34 39.66

A 60.34 39.66

A 60.59 39.41

A 60.59 39.41

The size of the 8 hours blend-shop problem when neglecting qualities (BS) is 670 continuous variables, 622 binary variables, 2,420 constraints with 370 as equality for the MILP. The conjugated NLP has 417 continuous variables, 491 constraints with 425 as equality. For the problem using factors (BSF), we have 754 continuous variables, 790 binary variables, 2,840 constraints with 401 as equality for the MILP. Table 7 shows the MILP factor results for each SS iteration in terms of the number of feasible solutions, objective (OBJ) function and the CPU resources. Table 7. Statistics of CPLEX and GUROBI solvers considering 8h horizon with 1h timestep in the BSF problem.


41


Page 42 of 56

CPLEX Number of solutions

OBJ (K$)

CPU (s)

GUROBI Number of solutions

OBJ (K$)

CPU (s)

st

1

354.20

0.19

3

354.20

0.67

nd

2

349.25

0.56

2

349.25

0.45

rd

2

348.00

0.70

3

348.00

0.47

th

2

347.59

0.20

2

347.59

0.46

th

2

347.44

0.25

2

347.44

0.77

1 SS 2 SS 3 SS 4 SS 5 SS

In Table 8, we scale the problem from a time-horizon of 8h (in Table 7) to 720h (30 days horizon). The solver GUROBI comparable to CPLEX is much faster in terms of CPU(s) with an average reduction of five times in the processing time. For the MILP problem using factors (BSF) in this case, we have 25,922 continuous variables, 24,480 binary variables, 73,436 constraints with 13,682 as equality (degrees-of-freedom equals to 50,402). The MILP-NLP gap in this 30 days horizon is the same presented in Table 5.

Table 8. Statistics of CPLEX and GUROBI for 720h with 1h time-step in the BSF problem. CPLEX Number of solutions

OBJ (K$)

st

39

nd rd th th

1 SS 2 SS 3 SS 4 SS 5 SS

CPU (s)

GUROBI Number of solutions

OBJ (K$)

CPU (s)

31,876.30

624.93

9

31,878.00

70.20

29

31,429.83

302.67

8

31,432.06

43.91

38

31,318.19

380.69

8

31,320.31

45.03

17

31,281.16

146.69

10

31,283.22

73.10

14

31,267.78

180.77

5

31,269.82

32.64

6. Industrial-sized case: gasoline blend scheduling The example 3 in Figure 11 illustrates an industrial-sized blend-shop for production of regular and premium grades of gasoline, BL-R and BL-P respectively. The properties to specify, treated


42

Page 43 of 56 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


as factors of quality, are octane number as ROAD and Reid vapor pressure (RVP). Although octane number calculation is highly nonlinear and depends on the influence of molecular interactions with intermediate streams known as synergistic and antagonistic effects, it is simply treated as a property that blends linearly or ideally with volume. In a closed-loop strategy, measurements for both inputs and outputs can be reconciled and estimated where updating flow and property biases in the mixture is a widespread practice to manage the slight differences between the real nonlinear blending and the simplified linear one.


43


Page 44 of 56

Figure 11. Example 3: industrial-sized blend scheduling of regular and premium gasoline.

For the MILP problem using factors for 1h as time-step and 180h as time-horizon, we have 17,247 continuous variables, 4,728 binary variables, 34,290 constraints with 9,411 as equality. In Table 9, we compare the CPU times of 1h and 2h time-step for CPLEX and GUROBI solutions with $ 628.215 as objective value.


44

Page 45 of 56 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


Table 9. CPLEX v. GUROBI for 7 days with 1h and 2h time-steps at 3% MILP relaxation gap. MILP Focus

CPU (seconds) CPLEX

GUROBI

1h

2h

1h

2h

feasibility

6539.8

360.5

923.6

95.6

optimality

18023.8

1420.8

460.2

42.5

Figure 12 shows the blended products in tanks PT1 to PT9 for the regular gasoline from BL-R and PT10 and PT11 for the premium gasoline from BL-P. In the blender modes of operations BL-R and BL-P can be verified in each time-window that the holdup or level of their respective tanks increases when there is an operation in the blender.

Figure 12. Example 3: industrial-sized blend scheduling of regular and premium gasoline considering 1h as time-step.


45


Page 46 of 56

7. Conclusion We have highlighted the application of blending with quantity, logic and quality details using a relatively simple and straightforward LP-based formulation by approximating nonlinear property pooling with constant but time-varying factors for qualities. The purpose of such optimization is to maximize profit and performance by increasing the sales of valuable products using the cheapest feedstocks subject to minimizing giveaway of key qualities explicitly in the MILP formulation albeit neglecting their nonlinear behavior which can be effectively polished or finetuned using the post nonlinear optimization. Finally, this technology has been implemented at an actual gasoline blend-shop where the factors are used to increase the likelihood of finding good qualogistics-feasible solutions in reasonable time with many benefit areas such as truly optimizing (and not just simulating) one- to two-months worth of blends into the future on-time and on-specification.

ASSOCIATED CONTENT *Supporting Information The configuration and flowsheet set of the illustrative case (Section 5.1) in the modeling platform IMPL (Industrial Modeling and Programming Language) is given. This material is available free of charge via the Internet at http://pubs.acs.org.

AUTHOR INFORMATION Corresponding Author *E-mail: [email protected].


46

Page 47 of 56 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


ACKNOWLEDGMENT The second author was funded by the Brazilian National Research Council (CNPq) under grant Ciências sem Fronteiras number 201542/2015-3.

NOMENCLATURE Subscripts e: units

,: unit-operations : in-port-states

66 : in-port-states downstream to out-port-states

: out-port-states

6 : out-port-states upstream to in-port-states

% : out-port-states of slacks or surpluses (factor-flow)

: properties, indices or qualities : time-periods Parameters ∆vNl :

∆vNwx :

minimum fill-draw delay maximum fill-draw delay

',, : factor (in volume flow) -,, : factor (in mass flow)

o : total number of time-periods

K̅, : lower bound for yield in the in-port-state of , at time K̅, : upper bound for yield in the in-port-state of , at time

K̅, : lower bound for yield in the out-port-state of , at time

K̅, : upper bound for yield in the out-port-state of , at time

bAg, : lower bound for utilization of the out-port-states of , at time bAg, : upper bound for utilization of the out-port-states of , at time


47


Page 48 of 56

: lower bound for utilization of the out-port-states of , at time bAg, bAg, : upper bound for utilization of the out-port-states of , at time

bj : lower bound for the uptime use , at time (for blender)

bj : upper bound for the uptime use , at time (for blender)

̅,, : volume-based property in the in-port-state of , at time

̅,, : volume-based specific gravity SG in the in-port-state of , at time

̅,, : lower bound for volume-based specific gravity SG in the in-port-state of , at time ̅,, : upper bound for volume-based specific gravity SG in the in-port-state of , at time ̅,, : lower bound for volume-based property in the out-port-state of , at time ̅,, : upper bound for volume-based property in the out-port-state of , at time

,, : mass-based property in the in-port-state of , at time

,, : lower bound for mass-based specific property in the out-port-state of , at time ,, : upper bound for mass-based specific property in the out-port-state of , at time ̅N, : lower bound for flow of the unit-operation , at time

̅N, : upper bound for flow of the unit-operation , at time

XXXN, ℎ : lower bound for flow of the unit-operation , at time Xℎ XXN, : upper bound for flow of the unit-operation , at time

̅,, : lower bound for flow of the stream between and at time ̅,, : upper bound for flow of the stream between and at time

Variables '& ,, : factor-flow (in volume) of the property p leaving the out-port-state % of , at time

-& ,, : factor-flow (in mass) of the property p leaving the out-port-state % of , at time ,, : unit-operation mass flow of , at time (No UOPSS)

,, : unit-operation mass flow of , at time (No UOPSS)

, : unit-operation volume flow of , at time (No UOPSS)

, : unit-operation volume flow of , at time (No UOPSS) N, : unit-operation flow of , at time


48

Page 49 of 56 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


ℎN, : unit-operation holdup of , at time

,, : unit-operation-port-state-unit-operation-port-state flow between and at time EN, : unit-operation setup (binary) of , at time

E,, : unit-operation-port-state-unit-operation-port-state setup (binary) between and at time

,, : volume-based property in the out-port-state of , at time

,, : volume-based specific gravity SG in the out-port-state of , at time ,, : volume-based property in the out-port-state of , at time ,, : mass-based property in the out-port-state of , at time

,, : mass-based property in the out-port-state of , at time dReN, : unit-operation startup of , at time

dRN, t: unit-operation switchover-to-itself of , at time

dRfN, : unit-operation shutdown of , at time REFERENCES

(1) Mouret, S.; Grossmann, I. E.; Pestiaux, P. A novel priority-slot based continuous-time formulation for crude-oil scheduling problems. Ind. Eng. Chem. Res., 2009, 48 (18), 8515-8528. (2) Castro, P.; Grossmann, I. E. Global optimal scheduling of crude oil blending operations with RTN continuous-time and multiparametric disaggregation. Ind. Eng. Chem. Res., 2014, 53 (39), 15127-15145. (3) Cerdá, J.; Pautasso, P. C.; Cafaro, D. C. Efficient approach for scheduling crude oil operations in marine-access refineries. Ind. Eng. Chem. Res., 2015, 54 (33), 8219-8238. (4) Menezes, B. C.; Kelly, J. D.; Grossmann, I. E. Phenomenological decomposition heuristic for process design synthesis of oil-refinery Units. Comput. Aided Chem. Eng., 2015, 37, 18771882.


49


Page 50 of 56

(5) Cerdá, J.; Pautasso, P. C.; Cafaro, D. C. A cost-effective model for the gasoline blend optimization problem. AIChE J., 2016, 62 (9), 3002-3019. (6) Cerdá, J.; Pautasso, P. C.; Cafaro, D. C. Optimizing gasoline recipes and blending operations using nonlinear blend models. Ind. Eng. Chem. Res., 2016, 55 (28), 7782-7800. (7) Kelly, J. D.; Menezes, B. C.; Engineer, F.; Grossmann, I. E. Crude-Oil Blend Scheduling Optimization of an Industrial-Sized Refinery: a Discrete-time Benchmark. In Foundations of Computer Aided Process Operations, 2017: Tucson, United States. (8) Haverly, C. A. Studies of the behavior of recursion for the pooling problem. Acm sigmap bulletin, 1978, 25, 19-28. (9) Mendez, C. A.; Grossmann, I. E.; Harjunkoski, I.; Kabore, P. A simultaneous optimization approach for off-line blending and scheduling of oil-refinery operations. Comput. Chem. Eng., 2006, 30, 614−634.

(10) McCormick, G. P. Computability of global solutions to factorable nonconvex programs: Part 1- convex underestimating problems. Math. Program., 1976, 10 (1), 147-175. (11) Al-Khayyal, F. A.; Falk, J. E. Jointly constrained biconvex programming. Math. of Oper. Res., 1983, 8 (2), 273-286. (12) Andrade, T.; Ribas, G.; Oliveira, F. A strategy based on convex relaxation for solving the oil refinery operations planning problem. In Foundations of Computer Aided Process Operations, 2017: Tucson, United States.


50

Page 51 of 56 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


(13) Assis, L. S.; Camponogara, E.; Zimberg, B.; Ferreira, E.; Grossmann, I. E. A piecewise McCormick relaxation-based strategy for scheduling operations in a crude oil terminal. Comput. Chem. Eng., 2017, 106 (2), 309-321. (14) Teles, J.P.; Castro, P.M.; Matos, H.A. Multiparametric disaggregation technique for global optimization of polynomial programming problems. J. Glob. Optim., 2011, 55, 227-251. (15) Baltean-Lugojan, R.; Misener, R. Piecewise parametric structure in the pooling problem: from sparse strongly-polynomial solutions to NP-hardness. J Glob Optim, 2017, 1-36. (16) Dey, S.; Gupte, A. Analysis of MILP techniques for the pooling problem. Oper Res, 63(2), 412-427. (17) Meyer, C.A.; Floudas, C.A. Global optimization of a combinatorially complex generalized pooling problem. AIChE J., 2006, 52(3), 1027–1037. (18) Misener, R., Floudas, C.A.: Advances for the pooling problem: modeling, global optimization, and computational studies. Appl Comput Math, 2009, 8(1), 3–22. (19) Misener, R., Floudas, C.A.: Global optimization of mixed-integer quadratically-constrained quadratic programs (MIQCQP) through piecewise-linear and edge-concave relaxations, Math Program, 2012, 136, 155–182. (20) Misener, R., Floudas, C.A.: ANTIGONE: algorithms for continuous integer global optimization of nonlinear equations, J Glob Optim, 2014, 59(2–3), 503–526.


51


Page 52 of 56

(21) Misener, R., Thompson, J.P., Floudas, C.A.: APOGEE: global optimization of standard, generalized, and extended pooling problems via linear and logarithmic partitioning schemes. Comput Chem Eng, 2011, 35(5), 876–892. (22) R. Karuppiah, I.E. Grossmann, Global Optimization for the Synthesis of Integrated Water Systems in Chemical Processes. Computers Aided Chemical Engineering, 2006, 30 (4), 650-673. (23) E. Ahmetovic, I.E. Grossmann, Global Optimization for the Synthesis of Integrated Water Systems in Chemical Processes. AICHE Journal, 2011,57(2), 434-457. (24) Kolodziej, S. P.; Grossmann, I. E.; Furman, K. C.; Sawaya, N. W. A discretization-based approach for the optimization of the multiperiod blend scheduling problem. Comput. Chem. Eng., 2013, 53, 122. (25) Harjunkoski, I., Martin, M. M. Modeling, Simulation and Optimization in the Process and Commodities Industries. In: Mariano M. Martin (Eds.) Introduction to Software for Chemical Engineers. (2015). CRC Press, Chapter 2, 11. (26) Castro, P. New MINLP formulation for the multiperiod pooling problem. AICHE J., 2015, 61, 3728-3738. (27) Castillo-Castillo, P. A., Mahalec, V. Inventory pinch gasoline blend scheduling algorithm combining discrete- and continuous-time models. Compt. Chem. Eng., (2016). 84, 611-626. (28) Lotero, I.; Trespalacios, F.; Grossmann, I. E.; Papageorgiou, D. J.; Cheon, M. S. An MILPMINLP decomposition method for the global optimization of a source based model of the multiperiod blending problem. Comput. Chem. Eng., 2016, 87 (6), 13-35.


52

Page 53 of 56 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


(29) Castillo-Castillo, P. A., Castro, P. M., Mahalec, V. Global Optimization of Nonlinear Blend-Scheduling Problems. Compt. Chem. Eng., (2017). 84, 611-626. (30) Chen, Y., Maravelias, C. T. (2017). Solution Methods for Multiperiod Blend Scheduling MINLP Models. In AICHE Annual Meeting, Minneapolis, MT: October. (31) Quesada, I.; Grossmann, I. E. Global optimization of bilinear process networks with multicomponent flows. Comput. Chem. Eng., 1995, 9 (12), 1219-1242. (32) Palacios-Gomez, F.; Lasdon, L.; Enquist, M. "Nonlinear Optimization by Successive Linear Programming", Manag. Sci., 1982, 28, 1106-1120. (33) Zhang, J.; Kim, N-H.; Ladson, L. An Improved Successive Linear Programming Algorithm, Manag. Sci., 1985, 31 (10), 1312-1331. (34) Kelly, J. D. The Unit-Operation-Stock Superstructure (UOSS) and the Quantity-LogicQuality Paradigm (QLQP) for Production Scheduling in The Process Industries. In Multidisciplinary International Scheduling Conference Proceedings: New York, United States, 327, 2005. (35) Zyngier, D., Kelly, J. D. UOPSS: A New Paradigm for Modeling Production Planning and Scheduling Systems. In European Symposium in Computer Aided Process Engineering: London, United Kingdom, 2012. (36) Menezes, B. C.; Kelly, J. D.; Engineer, F.; Grossmann, I. E.; Vazacopoulos, A. Generalized Capital Investment Planning of Oil-Refineries using MILP and Sequence-Dependent Setups. Comput. Chem. Eng., 2015, 80, 140-154.


53


Page 54 of 56

(37) Menezes, B. C.; Grossmann, I. E.; Kelly, J. D. Enterprise-Wide Optimization for Operations of Crude-Oil Refineries: Closing the Procurement and Scheduling Gap. Comput. Aided Chem. Eng., 2017, 40, 1249-1354. (38) Kelly, J. D.; Menezes, B. C.; Engineer, F.; Grossmann, I. E. Feedstock storage assignment in process industry quality problems. In Foundations of Computer Aided Process Operations, 2017: Tucson, United States. (39) Kelly, J. D.; Menezes, B. C.; Engineer, F.; Grossmann, I. E. Crude-oil blend scheduling optimization of an industrial-sized refinery: a discrete-time benchmark. In Foundations of Computer Aided Process Operations, 2017: Tucson, United States. (40) Kelly, J. D.; Menezes, B. C.; Grossmann, I. E. Decision Automation for Oil and Gas Well Startup Scheduling using MILP. Comput. Aided Chem. Eng., 2017, 40, 1399-1404. (41) Zyngier, D.; Lategan, J.; Furstenberg, L. A Process Systems Approach for Detailed Rail Planning and Scheduling Applications. Comput. Chem. Eng., 2018, 114, 273-280. (42) Kelly, J. D.; Menezes, B. C. Decision Automation for Lubes and Asphalt Production Scheduling using MILP with Sequence-Dependent Switchovers. In AICHE Spring Meeting, 2018: Orlando, United States. (43) Menezes, B. C.; Kelly, J. D.; Grossmann, I. E. Logistics Optimization for Dispositions and Depooling of Oil-Refineries: Closing the Production Scheduling and Distribution Gap. Comput. Aided Chem. Eng., 2018, 43, 1135-1140.


54

Page 55 of 56 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


(44) Franzoi, R. E.; Menezes, B. C.; Kelly, J. D.; Gut, J. W. Effective Scheduling of Complex Process-shops using Online Parameter Feedback in Crude-oil Refineries. Comput. Aided Chem. Eng., 2018, 44, 1279-1284. (45) Kelly, J. D., Zyngier, D. (2007). An Improved MILP Modeling of Sequence-Dependent Switchovers for Discrete-Time Scheduling Problems. Ind. Eng. Chem. Res., 46, 4964. (46) Zyngier, D., Kelly, J. D. (2009). Multi-product Inventory Logistics Modeling in The Process Industries. In: Wanpracha Chaovalitwongse, Kevin C. Furman, Panos M. Pardalos (Eds.) Optimization and logistics challenges in the enterprise. Springer optimization and its applications, 30, Part 1, 61.


55


Page 56 of 56

TOC GRAPHIC


56

Successive LP Approximation for Non-Convex Blending in MILP

Recommend Documents