《1.Introduction》

1.Introduction

The design and operation of process plants are subject to a number of uncertain parameters, such as fluctuating production outputs, ever-changing market prices and demands, and so on. A weakening product demand, catalyst deactivation, or financial woes for the corporation are but a few of the multitude of potential causes for optimizing the process. In order to ensure that the optimal process has been designed, it is imperative to certify that the proposed design operates both economically and safely, should physical or socioeconomic changes occur. Process synthesis [1] is a key aspect of plant design, having a significant financial impact on a firm. This has subsequently inspired substantial research in the field, in order to devise methods that will  accurately  design  not  only  the  optimal route, but also the most robust route. Furthermore, due to the computational complexity of the problem, it is essential that the formulated algorithms are efficient in order to minimize the computational cost. In general, the process synthesis problem for a chemical plant can be classified as an integrated system of the following major  components:  ① heat-exchanger  network  (HEN)  synthesis,  ② reactor path synthesis, and ③ separation systems synthesis. Reactor path synthesis is conducted to determine the optimal reactor configuration for a desired product. Within this problem, a number of factors ranging from the thermodynamic boundaries of the system to the reaction path are considered. There are three main techniques used for the optimization of reactor networks: the superstructure approach, the geometrical approach, and the combined targeting approach [2]. The superstructure approach requires a simultaneous synthesis of the entire process, which can be extremely effective, as all elements of the system can be determined with regards to the entire process. However, the validity of the solution is significantly dependent upon the detail of the model formulated to simulate the superstructure. The geometrical approach is founded upon the attainable region (AR) theory, which was developed to take into account the effect of having multiple interacting reactors and external heat exchangers [3]. The combined-targeting approach employs elements from both of the aforementioned techniques, as it stems from both the superstructure approach and the AR theory. Both the geometrical and the combined-targeting approaches incorporate the AR theory, limiting their utility. This is due to the fact that the AR theory integrates technical criteria, such as conversion and selectivity, as opposed to financial criteria; as a result, the solution may not be optimal from an economic standpoint. Bedenik et al. [2] have developed a novel technique to solve reactor network problems; this technique uses the superstructure approach and reformulates it to focus on the financial criteria of the model, including the impact of economic uncertainty over time, an aspect that will be discussed further in this report. Their technique is therefore more likely to be able to successfully find the optimal configuration, although it requires a thoroughly detailed mathematical model.

The vast majority of problems related to the field of process systems engineering involve continuous as well as discrete decisions. Typically, discrete decisions account for the underlying logic of the process, such as selection of a reaction pathway or expansion of capacity for planning problems, while continuous decisions account for product flow rates, amount of sales, and so forth. In order to mathematically model these kinds of decisions, continuous as well as 0–1 integer variables are employed, thus resulting in mixed-integer programming (MIP) problems. For a case in which the governing laws are expressed through nonlinear relations, the resulting problem is a mixed-integer nonlinear programming (MINLP) problem, which can generally be formulated as follows [1]:

where f is a scalar function; x is the nx-dimensional vector of continuous variables that belongs to the bounded set X; y is the ny-dimensional vector of binary variables that belongs to the discrete set Y; g is the vector of constraints that accounts for quality bounds, demand satisfaction, and so forth; and z* refers to the optimal solution of the optimization problem.

Because of their generic mathematical nature, MINLP problems have found application in a wide scope of problems such as process synthesis of reactor networks [4], HENs [5], process planning, and enterprise-wide optimization of manufacturing processes [6], to name just a few. Process synthesis forms a fundamental class of MINLP problems in the area of engineering, in which simultaneous decisions on the selection of processing units, interconnections, and design/operating variables need to be made [7]. Despite MINLPs’ ability to mimic the system under study more accurately than their linear counterparts (mixed-integer linear programming, MILP), these problems are difficult to solve because of the possible non-convexities that arise [8]. A number of numerical techniques have been proposed in the open literature to solve certain classes of MINLPs, such as the generalized Benders decomposition (GDB) [9], the outer approximation (OA) algorithm [10], and extended cutting plane (ECP) methods [11]. The aforementioned algorithms can be used in a rigorous manner under specific assumptions about the convexity of the objective function and the form of nonlinearities involved, while an iterative procedure is generally followed to typically converge to a (locally) optimal solution. Another class of MINLP solution algorithms falls under the scope of global optimization [12,13], in which specialized numerical techniques and convex approximations are employed so as to solve the corresponding problem within a tolerance of ε-optimality.

In principle, the solution of MINLPs can be computationally demanding, even  for  a  case  in  which  no  uncertainty  is  considered. However, mathematical models are susceptible to a number of uncertainties that can be broadly categorized as endogenous and exogenous [14]. Endogenous uncertainty is mostly encountered on the left-hand side of the constraints, such as reaction yields and stoichiometric coefficients, while exogenous uncertainty is usually located on the right-hand side (RHS) of the constraints and in the objective function coefficients (OFCs). In order to handle uncertainty within optimization problems, a number of solution  techniques  have been proposed in the literature; the main techniques are stochastic programming, robust optimization, and (multi-)parametric programming (mp-P) [15]. While the first two techniques require the availability of some form of data that can be used to characterize the uncertainty, such as probability distributions or type of uncertainty set, mp-P makes direct use of the mathematical model and, with optimization-based methodologies, aims  to  explicitly  characterize the effect of the uncertainty on the optimal solution. Through the solution of an mp-P problem, one generally aims to compute the optimization variables as explicit functions of the uncertain parameters, along with the corresponding regions—that is, critical regions (CRs)—of the parametric space where each function remains optimal. Fig. 1 provides a conceptual graph of such a problem.

《Fig. 1》

 

Fig. 1. Conceptual graph of the mp-P scope.

Having been studied for over 30 years [15,16], mp-P stands on solid theoretical foundations, with algorithms for every class of optimization problem.  However,  despite  the  constant  developments in the field of mp-P, the class of (multi-)parametric (mp)-MINLPs remains among the least studied, even for  a  case  in  which  only RHS uncertainty is considered. The reason for this lack of research is the non-convexity that is involved in the solution of the underlying parametric optimization problems. Early  works focused on the approximate solution of single parametric MINLPs (p-MINLPs)[17–21], while Dua and Pistikopoulos [22] studied the mp-MINLP case with convex objective functions and proposed a decomposition approach  in  which  successive  iteration  occurred  between  primal and master problems for the generation of valid upper and lower bounds, respectively. In their approach, the primal problem is an mp-nonlinear programming (mp-NLP) problem derived  by  fixing the related integer variables, while the master sub-problems aim to provide improved integer solutions. Dua et al. [23] studied the global optimization of mp-MINLPs. Polynomial parametric optimization techniques using algebraic geometry techniques were presented in Ref. [24], in which cylindrical algebraic decomposition and Gröbner bases theory are employed to solve the system of polynomial equations that arise in the model predictive control (MPC) scheme. Recently, Charitopoulos and Dua [25] presented a new algorithm for the exact solution of (multi-)parametric mixed-integer polynomial optimization (mp-MIPOPT) problems with emphasis on the explicit control of hybrid polynomial systems.

In the present work, we present two algorithms for the analytical solution of MINLPs for the case of logarithmic nonlinearities. The key idea is to analytically solve the square system of equations derived by the first-order Karush-Kuhn-Tucker (KKT) conditions using Gröbner bases theory and symbolic manipulation principles. The first algorithm is devised for the solution of deterministic MINLPs, while the second is an extension for the case of single-parametric perturbations on the RHS of the constraints of the problem. Problems related to process synthesis are then examined, following the proposed algorithms.

The remainder of the paper is organized as follows: In Section 2, the theoretical framework of the present work is established; next, two algorithms for the solution of a specific class of deterministic and parametric MINLPs are proposed. The proposed algorithms are then tested on two case studies related to process synthesis problems; these are presented in Section 3. After introducing the problems under study and their numerical optimal solution using stateof-the-art solvers, in Section 4.1, we validate the proposed algorithm for the deterministic case via a comparison with the solutions presented in Section 3. Next, in Section 4.2, we present the main computational steps as well as the results from the application of the proposed algorithm for the non-deterministic case. Finally, Section 5 contains concluding remarks.

《2.Theory and algorithms》

2.Theory and algorithms

Recently, Dua [26] presented an mp-inspired algorithm for the solution of MIPOPT problems. The algorithm is based on the analytical solution of the square system of equation derived by the first-order KKT conditions of the general MIPOPT problem. In general, a MIPOPT problem can be described as shown in Eqs. (2–6).

MIPOPT:

 

             

where x is the vector of continuous variables that belongs to the bounded set X; y represents the vector of integer variables that is ny-dimensional; h is the  vector of equality  constraints and  is nh-dimensional; g is the vector of inequality constraints and is ng-dimensional; and f is a scalar objective function. The integer variables of the problem (MIPOPT) are relaxed as they are continuous and are treated as parameters that are restricted within their respective bounds. Thus,  a  (multi-)parametric  polynomial  problem  (mp-PP)  arises,  as shown in Eqs. (7–11).

mp-PP:

 

 

As observed, the mp-PP is now parametric in the relaxed integer variables y. Eqs. (12–14) give the corresponding first-order KKT conditions.

KKT system:

where L ( x, y, μ, λ) = f ( x, y) + λT g ( x ) + μT h( x, y) is the Lagrange function of  the mp-PP and ng  denotes the number of  inequality constraints considered in problem mp-PP. Note that Eqs. (12–14) form a square system of polynomial equations that can be solved analytically with respect to the optimization variables and the Lagrange multipliers. Following the algorithm  that Dua [26]  proposed, the KKT system is solved analytically using symbolic manipulation software, thus computing all the possible solutions that satisfy Eqs. (12–14) as explicit functions of the y parameters. The explicit solutions are then validated through the evaluation of primal and dual feasibility conditions along with some constraint qualification, for all the possible combinations of the y vector.

In principle, the proposed methodology has the potential to be extended to more generic classes of MINLP problems for both the deterministic as well as the parametric case; this potential is examined next.

《2.1.Parametric programming-based algorithm for MINLPs》

2.1.Parametric programming-based algorithm for MINLPs

As mentioned in the previous section, Dua [26] studied the case of  MIPOPT  problems.  However,  since  the  key  component  of  the aforementioned  algorithm  is  the  analytical  solution  of  a  square system of equations through the use of Gröbner bases theory, the algorithm can be extended to classes of MINLP problems that involve logarithmic functions. By closely following the developments described by Dua [26], we devised Algorithm 1 for the case of logarithmic  nonlinearities.

《Algorithm 1》

Algorithm 1 Algorithm for the solution of MINLPs.

Note that in the present work, Steps 3–5 are performed in Mathematica 10 [27]. The linear independence constraint qualification (LICQ) is then evaluated upon obtaining the optimal solution, in order to confirm its optimality.

《2.2. MINLPs under uncertainty using parametric programming》

2.2. MINLPs under uncertainty using parametric programming

Algorithm 1 can facilitate deterministic MINLP problems. However, for a case in which uncertainty must be considered, the algorithm can be modified so as to handle parametric variabilities on the RHS of the constraints. The idea is to augment the vector of uncertain parameters; that is, in addition to the integer variables (y), a new vector of perturbations (θ) is considered. The resulting p-MINLP problem can be formulated as shown in Eqs. (15–20).

p-MINLP:

Note that θ is a bounded scalar parameter that is allowed to appear on the RHS of the equality/inequality constraints. Formulating the corresponding first-order KKT conditions results in Eqs. (21–23), which are now parametric in both y and θ.

p-KKT system:

Since both the parametric perturbation (θ) and integer variables (y) are considered as symbolic expressions, the parametric KKT (p-KKT) system remains a square system of nonlinear equations. The solution of the p-KKT system returns the candidate solutions, that is, the parametric explicit expressions of the optimization variables (x( y, θ)) and the Lagrange multipliers (λ( y, θ), μ( y, θ)). Note that even though the candidate solutions satisfy the square system of equations, they may violate the primal/dual feasibility conditions.

After the candidate solutions have been computed, the binary variables are fixed to their possible values; thus, the expressions of the optimization variables and the Lagrange multipliers are now explicit only in θ, that is, x(θ), λ(θ), and μ(θ). Note that in this step, the explicit expressions of the optimization variables are substituted back into the original constraints, that is, Eqs. (16, 17), and the primal feasibility of the candidate solutions is evaluated. Regarding the case of dual feasibility, we next evaluate the non-negativity condition of the Lagrange multipliers associated with the inequality constraints.

The candidate solutions that satisfy the primal/dual feasibility conditions are referred to as feasible. Each feasible solution has a parametric range of variability within which it remains feasible, that is, a CR. The main computational steps for the solution of p-MINLPs are given in Algorithm 2.

《Algorithm 2》

Algorithm 2 Algorithm for the solution of p-MINLPs.

Steps 1–3 are implemented in Mathematica 10, while Steps 4–5 as proposed by Vyas and Dua [28] are implemented in Excel. The generated output from Mathematica 10 should create objective functions in terms of the binary variables and the uncertain parameters. These can then be used to help determine the optimal configuration for the system, by clearly showing how varying parameters will impact the objective function.

《3.Case studies》

3.Case studies

This section presents two case studies related to process synthesis to illustrate the proposed methodology. First, a short description of the case study is given, along with the related mathematical model, which results in a MINLP problem with logarithmic nonlinearities. Next, we report the numerical solution of each problem using a state-of-the-art numerical optimization solver. For the solution of the corresponding MINLPs, the SBB solver is employed, with CONOPT3 as the NLP sub-solver. The problems are implemented in GAMS 24.4.1 in a Dell workstation with a 3.7 GHz processor, 16 GB RAM, and a Windows 7 64-bit operating system.

《3.1. Case Study 1》

3.1. Case Study 1

Case Study 1 is adapted from the book by Floudas [1]. This process synthesis problem relates to the production of component C from raw materials A and B via Processes 1, 2, and 3. Product C can only be produced through one of three routes: through Process 1 only, through Processes 1 and 2, or through Processes 1 and 3. Processes 2 and 3 cannot take place simultaneously. Fig. 2 gives a conceptual representation of the related process flowsheet.

《Fig. 2》

Fig. 2. A process flow diagram showing the flowsheet to be optimized.

The objective is to minimize the cost minus the revenue. The corresponding optimization problem is formulated as a MINLP involving Eqs. (24–38).

Table 1 gives the optimal solution, as reported in Ref. [1].

《Table 1》

Table 1 Optimal solution for the example from Ref. [1].

《3.2. Case Study 2》

3.2. Case Study 2

The second case study is adapted from Ref. [29], and is a modified version of Example 1 found in that paper. The formulation of the problem is as follows:

where U = 2 is a “big-M” type of parameter. Table 2 gives the optimal solution, as reported in Ref. [29].

《Table 2》

Table 2 Optimal solution for Case Study 2.

This formulation is modified by relaxing the upper bound to the variables,  thereby  ignoring  Eqs.  (48–50);  the  optimal  solution  is then computed using GAMS 24.4.1. The optimal solution for this example is obtained using the SBB solver, as shown in Table 3.

《Table 3》

Table 3. Optimal solution for modified Case Study 2.

《4.Results》

4.Results

The previous section illustrated the process synthesis problems  under study and reported their optimal numerical solutions. This section now illustrates the main computational steps of the proposed algorithms. First, in Section 4.1, the deterministic case of each problem is solved and validated with the corresponding numerical solution, as computed in the previous section. Next, in Section 4.2, the parametric optimization problems are introduced and solved, employing the algorithm for p-MINLPs that is proposed in this study.

《4.1. Parametric programming-based algorithm for MINLPs》

4.1. Parametric programming-based algorithm for MINLPs

4.1.1. Case Study 1

After the binary variables in the model have been considered as symbolic parameters, the first step of the algorithm is to formulate the Lagrange function and the first-order KKT equations. The Lagrange function is given by Eq. (52) and the corresponding first-order KKT conditions are given by Eqs. (53–70).

Eqs. (53–70) are subsequently solved using Mathematica 10, with the binary variables being considered as parameters. The output of this step is a list of candidate solutions for which the optimization variables and the corresponding Lagrange multipliers are given as explicit functions of the binary variables. In total, 98 candidate solutions are computed. However, 40 of the solutions are ignored, as they involve contradicting requirements in terms of primal/dual feasibility; the remaining 58 are investigated further. Table 4 presents a sample of the output generated by Mathematica 10.

《Table 4》

Table 4 Sample of the candidate solutions computed after Step 3 of Algorithm 1 for Case Study 1.

As shown in Table 4, the optimization variables and the Lagrange multipliers are, in general, explicit expressions of  the binary variables. Also note that, through a preliminary screening test with respect to the dual feasibility of the candidate solutions, Solutions 1, 3, and 4 can be removed from further consideration, as they violate  the  requirement  for  non-negative  Lagrange  multipliers.  Next, following Algorithm 1, all the binary variables are fixed according to the possible combinations. There are three binary variables, resulting in seven possible combinations for the process, as shown in Table 5.

《Table 5》

Table 5 Possible combinations of binary variables for Case Study 1.

Each solution is evaluated by fixing the binary variables to the combinations shown in Table 5. This step can be coded and implemented in Mathematica 10 using the built-in command “Reduce.”

On examination of Table 6, the vast majority of the candidate solutions computed by Mathematica 10 are found to be infeasible, as only integer combinations of  four candidate solutions remain, with ranging values of the objective function. It is important to note that even if a solution is feasible, it may not be a viable configuration. This is due to the fact that in order for the process to be valid, it must be able to produce the desired product. Evidence of this fact can be seen when looking at Solution 12, which has the integer solution of y1 = 0, y2 = 1, and y3 = 1; there are no flow rates, and although y1  is the unit that is meant to generate the final product, it is not operating in this instance.

《Table 6》

Table 6 Final integer feasible solutions for Case Study 1.

Table 6 shows that the global minimum is –1.92. This value occurs when [y1, y2, y3] = [1, 0, 1]; it is derived from two different symmetrical candidate solutions. As expected, the solution computed by the proposed algorithm is identical to the one reported in the book by Floudas [1], thus validating the solution. Finally, in order to evaluate with the LICQ, it is important to determine which constraints from the original model are active at the optimal point. At the optimal point, the flow rates of B2  and A2  are zero, meaning that Eqs. (25, 29) are not considered. Therefore, the constraints considered for the LICQ are Eqs. (26–28, 30). The nonlinear equation, Eq. (26), must be linearized at the optimal point. Using the Taylor approximation at the optimal point yields:

where * refers to the value of the variable at the point of linearization. A matrix of the system of equations is then formulated:

Whether the system adheres to the LICQ can now be ascertained by  calculating  the  determinant  of  the  square  matrix.  When  the determinant of a matrix is zero, the equations in the system are dependent on each other, suggesting that there could be multiple solutions. However, if the determinant is non-zero, the lines are independent and will only intersect once, at a single solution. In this example, the determinant is 0.429. This means that the system is linearly independent, and the optimal solution has been determined.

4.1.2.Case Study 2

For the second case study, the Lagrange function and first-order KKT conditions are given by the following equations:

Mathematica 10 is then employed to solve Eqs. (74–84), in which the binary variables are considered as parameters. The output generates a total of 59 candidate solutions; however, similar to Case Study 1, 17 of them are removed from further consideration because of violations of the dual feasibility conditions, thus resulting in 42 possible solutions. Table 7 presents a sample of the output generated from Mathematica 10.

《Table 7》

Table 7 Sample of the candidate solutions computed after Step 3 of Algorithm 1 for Case Study 2.

Like Case Study 1, this case study has three binary variables; thus, the output from Mathematica 10 is subjected to the same variation of binary variables, as shown in Table 5. Table 8 gives the final integer feasible solutions.

《Table 8》

Table 8 Final integer feasible solutions for Case Study 2.

On examination of Table 8, the global minimum is seen to be 5.58 with an integer vector when [y1, y2, y3] = [0, 1, 0]. The same solution was computed from SBB/CONOPT3, as shown in Table 3. This result reinforces the applicability of the algorithm to MINLP problems. Finally, Table 9 presents a computational comparison of the proposed algorithm with the commercial solvers DICOPT, SBB, BARON 16.3, and ANTIGONE 1.1. As seen in Table 9, the proposed algorithm is out- performed for the deterministic instance by the rest of the solvers, with DICOPT being the solver that converges the fastest. However, note that for the non-deterministic case, a series of MINLP problems would have to be solved in order to obtain the solution of the corresponding p-MINLP. In contrast, when following the proposed methodology, the computational time remains the same as with the deterministic case.

《Table 9》

Table 9 Computational comparison of the performance of the proposed algorithm for the solution of the deterministic instance of Case Study 2.

《4.2.MINLPs under uncertainty using parametric programming 》

4.2.MINLPs under uncertainty using parametric programming 

4.2.1.Case Study 1

The uncertain variables to be investigated in this example are B1 and C. This is due to the fact that B1 is the flow rate of fresh feed brought into the system, which will subsequently be converted into the desired product, and C is the product of the system, which carries the possibility of varying due to consumer demand.

4.2.1.1. B1  as the uncertain variable

To investigate the effect of B1 on the process, Mathematica 10 is used to solve the following first-order KKT conditions, where the binary variables and B1  are considered as parameters.

The first-order KKT conditions shown above have been altered from the formulation derived in Eqs. (53–70), which was used when optimizing the process in the previous section. The non-negativity constraint from the original model and the Lagrange derivative of B1 are not included in this formulation, since B1  is considered to be a parameter. This results in an output that consists of explicit solutions in terms of the binary variables and B1. Solving the system of Eqs. (86–101) results in 39 candidate sets of explicit solutions; however, 13 of these violate dual feasibility conditions and, as such, are removed from further consideration.

The greatest and optimal exit flow rate of C from the process is 1, as can be seen in Eq. (28) and Table 1, respectively. Assuming that the reaction from B to C is 1:1, the flow rate of B1 is varied from 0 to 1 in order to examine its impact upon the process. The first step of the analysis is to assess the feasibility of each solution. This was achieved using Excel by varying the flow rate of B1 from 0 to 1 and by varying the configuration of binary variables as shown in Table 5. This was done to ensure that the flow rates and Lagrange multipliers of the inequality constraints are positive.

The objective function at each of the solutions shown in Table 9 must then be examined, as B1 varies. Table 10 provides the objective function for each solution, as derived by Mathematica 10.

《Table 10》

Table 10 Optimal explicit results when B1 is uncertain for Case Study 1.

The functions shown in Table 10 are then evaluated at their feasible configurations, as shown in Table 11, through the corresponding integer combinations, while manipulating B1 to examine its impact on the optimal solution. Table 12 provides the results.

《Table 11》

Table 11 Binary combinations when B1 is the uncertain parameter in Case Study 1.

《Table 12》

Table 12 Objective value for each integer feasible solution under varying values of B1  for Case Study  1.

An examination of Table 12 clearly shows that the lowest value of the objective function is –1.92 when B1 is equal to 0 in Solution 1, with a configuration of y1 = 1, y2 = 0, and y3 = 1. This result yields the same optimal solution as can be seen in the literature (Table 1). More  importantly,  however,  as  the  flow  rate  of  B1 increases,  the configuration with the lowest value of the objective function changes. This result strongly suggests that changing variables can have a significant impact on the optimal process route within a plant. The same result can be seen more clearly in Fig. 3, which shows how the optimal configuration changes as the flow rate of B1  increases.

《Fig. 3》

 

Fig. 3. Explicit optimal objective function for the range of variability of the uncertainty parameter B1.

To verify that the optimal configuration of the plant changes depending on an uncertain parameter, the deterministic formulation was solved again, with B1 as a varying parameter. Table 13 provides the corresponding results.

《Table 13》

Table 13 Effect of variations in the value of B1 on the optimal solution for Case Study 1.

Table 12 underlines the fact that the optimal configuration of a process can change based on uncertain parameters. Both techniques—that is, the technique using Mathematica 10 and the screening process in Excel and the numerical solution of the MINLP problems—generated the same output. This result reiterates the validity of the algorithm used for the process synthesis. It confirms the suggestion that the optimal flow rate of B1 is 0, with a configuration of y1 = 1, y2 = 0, and y3 = 1. It also demonstrates that an increase in the flow rate of B1 leads to the process becoming less profitable. In turn, this result implies that, should the need arise for the purchase of B1 from the market, the process will eventually operate at a loss and will not be financially viable.

4.2.1.2.C as the uncertain variable

Next, we investigate the effect of C on the process. Algorithm 2 is  again  employed,  and  the  corresponding  p-KKT  conditions  are formulated  and  solved  using  Mathematica  10.  Solving  the  p-KKT results in 32 candidate solutions, with 11 solutions removed during preprocessing; thus, 21 candidate solutions are considered for the subsequent algorithmic steps. Due to the constraint in Eq. (28), the maximum possible flow rate of C is 1. Therefore, to examine how C impacts the process, this variable will be varied from 0 to 1.

The objective function at each of the solutions shown in Table 13 must then be examined, as C varies. This can be achieved by examining the functions, which are derived by Mathematica 10 for each solution, of the objective function in terms of the binary variables and the uncertain variable, as shown in Table 14.

The functions shown in Table 14 can then be evaluated using the binary combinations shown in Table 15 while manipulating C; the values of the objective function are reported in Table 16.

《Table 14》

Table 14 Optimal explicit results when C is uncertain for Case Study 1.

《Table 15》

Table 15 Binary combinations when C is the uncertain parameter in Case Study 1.

《Table 16》

Table 16 Objective value for each integer feasible solution under varying values of C for Case Study 1.

An examination of Table 16 clearly shows that the lowest value of the objective function is –1.92 when C is equal to 1, with the integer solution of y1 = 1, y2 = 0, and y3 = 1. This binary combination yields the same optimal solution that was reported in the literature, as shown in Table 1. This result reinforces the solution that is obtained when the process is optimized parametrically and under uncertainty, and when B1 is the uncertain parameter. It also reiterates the fact that an uncertain variable can influence the optimal process route. Fig. 4 shows how the optimal configuration changes as the flow rate of C increases.

《Fig. 4》

Fig. 4. Explicit optimal objective function for the range of variability of the uncertainty parameter C for Case Study 1.

To verify that the optimal configuration of the plant changes depending on the uncertain parameter, the original formulation—Eqs. (24–38), where Eq. (31) is removed to ensure that C is considered as a parameter—was entered into GAMS 24.4.1.

4.2.2. Case Study 2

The uncertain variable to be investigated in this example is x1. Mathematica 10 is used to solve the following p-KKT conditions, in which the binary variables and x1  are considered as parameters.

Using Mathematica 10, 23 candidate solutions were computed, with eight solutions being removed during preprocessing. The remaining 15 solutions were analyzed further. Based on Eq. (48) from the original formulation, the maximum value for x1 was 2. Therefore, to study how x1 affects the process, this variable will be varied from 0 to 2, and Excel will be used to determine whether each solution is feasible at each configuration of binary variables.

It is important to emphasize that the feasibility of the solution can change, depending on the function of each of the variables. Each variable in the solution may be impacted differently by the value of the uncertain parameter. This can subsequently lead λ to fluctuate in a given range, potentially leading to infeasibility. The objective function at each of the integer solutions shown in Table 17 must then be examined as x1 varies. This can be achieved by examining the objective functions derived by Mathematica 10 for each solution, which are provided in Table 18.

The functions shown in Table 18 can then be evaluated at their feasible configurations, as shown in Table 17, while manipulating x1 to examine its impact on the optimal solution. Table 19 presents a sample of the impact that variations of x1  have on the optimal explicit solution.

《Table 17》

Table 17 Integer feasible solutions and parametric ranges for varying x1 for Case Study 2.

《Table 18》

Table 18 Optimal explicit results when x1 is uncertain for Case Study 2.

《Table 19》

Table 19 Different solution instances for two integer feasible solutions under varying  x1 for Case Study 2.

From Table 19, it is evident that the lowest value of the objective function is 5.583 when x1 is equal to 1.8 from Solution 2, with a configuration of y1 = 0, y2 = 1, and y3 = 0. This result gives an extremely similar optimal solution to the one provided in Table 3. Upon increasing the accuracy of the variation of x1 to two decimal places, the optimal solution reported in Table 3 is achieved. This result reiterates the fact that the algorithm can be applied to p-MINLP problems in order to obtain optimal solutions. It also reinforces the fact that an uncertain variable can influence the optimal process route. This fact is demonstrated more clearly in Fig. 5, which shows how the optimal configuration changes as the value of x1  increases.

To verify the accuracy of Fig. 5, the original formulation, Eqs. (39–51), was solved using GAMS 24.4.1 in order to ensure that x1 is considered as a parameter.

《Fig. 5》

Fig. 5. Explicit optimal objective function for the range of variability of the uncertainty parameter x1 for Case Study 2.

《5.Conclusions》

5.Conclusions

Under the constant pressure of uncertainty, the need for efficient and reliable decision-making in the process industry is increasingly important. The problem of decision-making in the process industry is indubitably a multi-layered issue, and decisions made at an early stage, such as those related to process synthesis, can determine the sustainability of the process. Current developments in the design of new-generation bio-refineries [30,31] advocate for the importance of uncertainty-aware decision-making, especially for the stage of process synthesis.

Motivated by this situation, in this study, we introduce two algorithms for solving deterministic and parametric MINLPs, respectively, for the special case in which the nonlinearities are expressed by logarithmic functions. The idea of the proposed algorithm is two-fold: First, the binary variables are treated as uncertain parameters and restricted within their respective bounds; next, the square system of equations derived from the first-order KKT conditions is solved analytically using a symbolic manipulation technique. This process results in a set of candidate solutions that are examined further so as to only keep the optimal solution. For the case of single-parametric uncertainty,  a  modified  version  of  the  aforementioned  algorithm is proposed, in which the varying parameter and binary variables are treated together as uncertain parameters. Following the proposed parametric algorithm, the optimizers—and consequently the objective function—are computed as exact explicit functions of the uncertain parameter together with the corresponding parametric ranges, for which each explicit expression remains optimal. Despite the merits of the proposed algorithm, there is a trade-off between the ability of the algorithm to compute the exact explicit solution and the size of the problems handled. Current developments in our group aim toward a computationally more efficient implementation of the proposed algorithms.

《Acknowledgements》

Acknowledgements

The authors gratefully acknowledge financial support from EPSRC grants (EP/M027856/1, EP/M028240/1).

《Compliance with ethics guidelines》

Compliance with ethics guidelines

Vassilis M. Charitopoulos, Lazaros G. Papageorgiou, and Vivek Dua declare that they have no conflict of interest or financial conflicts to disclose.