Decentralized control using selectors for optimal steady-state operation with changing active constraints

We study the optimal steady-state operation of processes where the active constraints change. The aim of this work is to eliminate or reduce the need for a real-time optimization layer, moving the optimization into the control layer by switching between appropriately selected controlled variables (CVs) in a simple way. The challenge is that the best CVs, or more precisely the reduced cost gradients associated with the unconstrained degrees of freedom, change with the active constraints. This work proposes a framework based on decentralized control that operates optimally in all active constraint regions, with region switching mediated by selectors. A key point is that the nullspace associated with the unconstrained cost gradient needs to be selected in accordance with the constraint directions so that selectors can be used. A main benefit is that the number of SISO controllers that need to be designed is only equal to the number of process inputs plus constraints. The main assumptions are that the unconstrained cost gradient is available online and that the number of constraints does not exceed the number of process inputs. The optimality and ease of implementation are illustrated in a simulated toy example with linear constraints and a quadratic cost function. In addition, the proposed framework is successfully applied to the nonlinear Williams–Otto reactor case study.


Introduction
The integration of optimization and control is very important when designing the control system for a process.The main objective of the control system is to keep the process stable and operating at the economically optimal operating point.Although these two objectives can be assessed simultaneously, for example, using economic model predictive control (EMPC) [1], a simpler, and in most cases equally optimal,1 approach is to decompose the system hierarchically into an optimization and a control layer as shown in Fig. 1, where setpoints   are used to connect the two layers.The setpoints may need to be updated due to disturbances that affect the process economics.In the standard implementation in Fig. 1, the real-time optimization (RTO) and setpoints update is performed on a slow time scale based on a detailed nonlinear process model and the estimated states of the process.In most cases, the RTO layer is static.
Based on the concept of Morari et al. [2] of feedback optimizing control, the aim of the current paper is to move the real-time optimization, or at least parts of it, into the control layer.A recent review on this topic is given in Krishnamoorthy and Skogestad [3], where the authors state some of the challenges with RTO implementation, including the cost of developing the model, the uncertainty related to the model and its parameters (or disturbances), and human aspects related to the maintenance of an optimization layer in addition to the already existing digital control system (DCS).The importance of feedback optimizing control lies in being able to reject disturbances that affect economic performance in a simple manner, without relying on an upper optimization layer that may sometimes not even exist.To that end, an appropriate selection of the controlled variables (CVs) for the control layer is important.This is the main idea of self-optimizing control [4].It is particularly important to include the active constraints as CVs, that is, the constraints that are optimally at their limiting value [2,5].If information about the cost gradient is available, the optimal CVs are the active constraints plus the reduced cost gradients, and by controlling these at a constant setpoint of zero we may eliminate the optimization layer [6].This choice of CVs is valid if the set of active constraints does not change in the considered operating region.https://doi.org/10.1016/j.jprocont.2024.103194Received 5 October 2023; Received in revised form 24 February 2024; Accepted 3 March 2024 Fig. 1.Standard optimizing control implementation with separate layers for realtime optimization (RTO) and control (, which can be ..MPC or PID). denotes the (economic) cost function to be minimized,  the process model,  the process constraints,  the model states,  the disturbances, and  the process inputs (MVs).
Dealing with changes in active constraints has been a concern in previous works.For example, Cao [7] implemented a cascade control structure with selectors to avoid constraint violation by the lower self-optimizing layer, and Graciano et al. [8] applied MPC with zone control to the same end.A global self-optimizing control method for changing active constraints has been proposed by Ye et al. [9], where the goal is to minimize the average loss obtained with a single set of CVs.However, in a new active constraint region, not only do the active constraints change, but the directions related to the reduced cost gradient change accordingly.This means that to eliminate the RTO layer one needs to change the control layer in Fig. 1 during operation, both in terms of the selected CVs and the corresponding feedback controller .With this perspective, Manum and Skogestad [10] has considered a centralized, steady-state analysis on switching control structures, with different CVs for each region.
However, the implementation of such a region-based control strategy quickly becomes impractical.This is because the number of active constraint regions grows exponentially with the number of constraints.Let   denote the number of process inputs or manipulated variables (MVs) and   the number of independent constraints.The upper bound on the number of active constraint regions is 2   , which is reached when all constraint combinations are feasible [11].In each region, we ideally need a new controller , and if we want to use decentralized control then we need to design   single-input single-output (SISO) controllers in each region.For example, with   = 4 and   = 5, there could be up to 2 4 = 16 constraint regions, which may require the tuning of 2   ⋅   = 16 ⋅ 5 = 80 SISO controllers.Even though some CVs are reused between regions, the number of necessary SISO loops will be high.
The key contribution of this paper is to propose a simple and generic region-based control structure with only   +   SISO controllers, as represented in Fig. 2, with the same set of unconstrained variables ( 0 and  0 ) in all operating regions. 2 Considering the previous example, this structure would have only 4 + 5 = 9 SISO controllers.In the paper, we show that the unconstrained variables are obtained from   projections of the full cost gradient with respect to the inputs, ∇   .This leads to   gradient controllers ( 0 and  0 ) and   constraint controllers (  ).However, at any given time, only a subset with   of the   +   controller outputs is implemented as process inputs, with the switching logic choosing between the controller outputs  0 and   .Fig. 2. Proposed optimizing control implementation, assuming   ≥   .The controllers  0 ,  0 and   are usually single-variable PID controllers.The projection (nullspace) matrices  0 and  are defined in Eq. ( 7) and Eq. ( 8), respectively.There is no CV 0 ,  0 , and  0 if   =   .Note that the optimization layer in Fig. 1 is eliminated, and an estimate ∇  Ĵ of the cost gradient is needed.The switching logic takes care of the change between active constraint regions.In this paper, this logic is decentralized to   individual blocks, see Fig. 3, which can be implemented as min or max selectors according to Theorem 3.
The second key contribution of this paper is to show that the switching logic in Fig. 2 can be effectively implemented using   min or max selectors, which are well-known advanced control elements and commonly used in practical control applications.An important decision is to pair each constraint to an MV, but this pairing problem is not addressed in this paper (the interested reader is referred to Skogestad and Postlethwaite [12]).The main assumptions in this work are that we have at least as many MVs as constraints (  ≥   ), and that an estimator for the unconstrained cost gradient ∇   is available.In terms of cost gradient estimation, there are several methods available (see Krishnamoorthy and Skogestad [3]), and in this work, we use the simple model-based approach of dynamic state estimation and model linearization proposed by Krishnamoorthy et al. [13].
Selectors have been used in industry to switch between CVs since the 1940s [14].Selectors are also used in academic case studies on optimal operation [11,15].In these case studies, a control structure is proposed for the nominal operating region, with added logic elements and control loops to deal with the neighboring regions.However, the treatment of the unconstrained degrees of freedom is not clear.Krishnamoorthy and Skogestad [16] proposes a framework for constraint handling using min and max selectors, focusing on systems with a single MV, and therefore not considering the changes of reduced gradients for the unconstrained variables.To the best of the authors' knowledge, even though a general scheme for the paradigm of region-based control is proposed in the review paper by Krishnamoorthy and Skogestad [3], a systematic procedure for designing a decentralized control structure for optimal operation of generic multivariable systems has not yet been explored, as well as whether there are any fundamental limitations for the design of such systems.In this work, we explore these topics, and we describe a class of multivariable systems for which a decentralized control structure is always possible.

Decentralized control framework for optimal operation
We consider a generic, steady-state optimization problem given by: is the function that returns the vector of inequality constraints,  ∈ R   is the vector of decision variables (MVs), and  ∈ R   is the vector of disturbances.Note that the states  (see Fig. 1) have been formally eliminated from the equations, such that  and  are functions only of the independent variables  and .Introduce the Lagrange function (, , ) =  (, ) +   (, ).Then, for a given value of , define  * as the solution of Eq. ( 1), which satisfies the Karush-Kuhn-Tucker (KKT) conditions [17]: Here,  is the vector of Lagrange multipliers associated with the inequality constraints, and  * is its optimal value.We remark that the KKT conditions only imply that the solution is a stationary point, and they are also satisfied by local minima or maximum and saddle points.We do not address these issues in this work, and we consider that the optimization problem in Eq. ( 1) is convex.While these optimization problems can be efficiently solved using numerical methods, we here focus on how to solve these problems with feedback control.For this, we rewrite the KKT conditions as control objectives, which allows us to embed the optimization into the control layer design.
The set of active constraints  is defined as the set that satisfies   ( * , ) = 0 for  ∈ .For convenience, define   ∶ R   × R   → R   as the function that returns the active constraints.Define the matrix: as the gradient of the constraints with respect to the MVs, and the matrix    = ∇    (, ) as the gradient of the active constraints with respect to the MVs.If the set of active constraints  is known, Jäschke and Skogestad [6] prove that optimality can be attained by controlling to zero the active constraints and the associated reduced cost gradient.Their result is given by the following theorem: Theorem 1 (Optimal Controlled Variables).Consider the optimization problem in Eq. (1), where we assume that linear independence constraint qualification (LICQ) holds.We assume that the set of optimally active constraints  is known.Let   ∈ R   ×(  −  ) be a basis for the nullspace of    such that: Further, define the reduced cost gradient as: Then controlling   (, ) = 0 and ∇ ,  (, ) = 0 results in optimal steady-state operation.
Proof ([6]).If the active constraints  are known, the necessary optimality conditions (2) are equivalent to: where  *  > 0 is the optimal vector of Lagrange multipliers for the active constraints.Premultiplying ∇  ( * , ) by    leads to: Since by definition      = 0, the optimality conditions are equivalent to   ( * , ) = 0 and    ∇   ( * , ) = 0, which are   equations that fully determine  * because   is full rank, and the associated optimal Lagrange multiplier can always be found as ) −1    ∇   ( * , ).Therefore, enforcing   ( * , ) = 0 and    ∇   ( * , ) = 0 leads to satisfying (6), which is equivalent to satisfying (2).□ In terms of feedback control, Theorem 1 says that   and ∇ ,  (both with setpoints 0) are the steady-state optimal CVs for a given operating region where the active constraints do not change.Here, the reduced cost gradient ∇ ,  =    ∇   is defined as the gradient in the unconstrained directions as given by the nullspace   of the active constraints [6].If the system is to operate at another active constraint region, however, the CVs need to change, and if shifts in operating regions happen in real-time, the control system needs to automatically detect these region switches.The main idea of this work is to design a decentralized control structure, see Fig. 2, for all possible active constraint regions of the optimization problem in Eq. ( 1).The main assumption for guaranteeing the existence of this decentralized control structure is as follows: The matrix   is always full row rank, and the number of constraints is not greater than the number of MVs, that is, rank(  ) =   , and   ≥   .
This not only guarantees LICQ for any set of constraints that may be optimally active, but it also guarantees the existence of decoupled CVs for optimal operation, as shown in the next theorem.For use in the next theorem, define  0 as an orthonormal basis of the nullspace of   , that is: The matrix  0 represents the unconstrained directions that are never in conflict with constraint control.Note here that  0 is an empty matrix (nonexistent) if we have as many constraints as inputs (  =   ).Further define   − as the matrix containing all but the th row of   , and define: as a matrix of   columns, where each column   is a unitary vector such that: Each vector   represents the direction that may conflict with the corresponding constraint   , as shown next.
Theorem 2 (Optimal Switching Between CVs).Given that Assumption 1 holds and that the active constraint index set is , the following control strategy allows for optimal operation: Proof.To prove Theorem 2, it is sufficient to prove that the controlled variables are equivalent to the necessary first-order optimality conditions.Firstly, it is useful to note that, due to its construction,     = (     ) ê , with    being the th row of   , and ê being the th unit vector from the standard basis.Additionally, if the active constraint set is , and the inactive constraint set is  = {1, … ,   }−, the optimality conditions can be written as: Let   be the matrix with columns equal to   for  ∈ .Then, premultiplying ∇   by [    0 ]  leads to: Here, (     )   * = 0, because     = (     ) ê , and, from the optimality conditions, ê   * =  *  = 0 for  ∈ .Therefore, the optimality conditions become which are the CVs proposed in addition to   (, ) = 0 for  ∈ .Similarly to Theorem 1, this fully defines the operational degrees of freedom, and a suitable vector of Lagrange multipliers can be found.□ Note that the matrix used in the proof of Theorem 2 is a particular parametrization of the nullspace matrix   from Theorem 1, and therefore both results are equivalent for a given active constraint region.In Theorem 2 however, we specify an ideal association between CVs such that the handling of region switching may be done in a decentralized fashion, avoiding changes in the rest of the control structure.For instance, if the th constraint changes from inactive to active, only the corresponding unconstrained degree of freedom  0  =    ∇   will become uncontrolled, and the remaining CVs are kept unaltered.In addition, the matrix is designed to be full row rank, and therefore all operational degrees of freedom are filled for any active set .The choice of building the vectors   unitary and orthogonal to  0 is purely for the uniqueness of the solution, as one could propose another projection  ′  =   +  0  for any nonzero scaling factor  and any vector , and optimal operation would still be attained, as  0  ∇   is always optimally zero.
Theorem 2 states a general set of feedback control objectives to attain optimal operation.It does not specify the type of controller to be used, and one may apply these results to obtain optimal operation with conventional tracking MPC with switching objectives to eliminate the RTO layer.This would be useful for cases where decentralized control performs poorly, but one still wishes to propose a simple control layer.In this work, however, we choose to explore the implications of this result for decentralized control, which is often more easily implemented in practice.
Pairing of MVs and CVs.It should be noted that Theorem 2 makes no distinction about the pairing between MVs and CVs, and it is left for the practitioner to make this pairing taking into account controllability and performance aspects.However, the theorem states the optimal association between CVs for region switching, which means that the control of    =   and  0  =    ∇   must be performed by the same MV in the case of a decentralized framework.From now on, it is considered that the MVs are ordered such that   is used to control the pair    =   and  0

𝑖
for  ≤   .These considerations lead to the control structure presented in Fig. 3. Here, ∇  Ĵ represents the estimate of the cost gradient (∇   ),    represent the individual constraint controllers,  0  represent the individual gradient controllers that are conditionally active (..only one of    and  0  is active at any given time), and  0  represent the individual gradient controllers that are always active.It is important that the controllers    and  0  include anti-windup action so that the integral modes in the inactive controllers do not grow indefinitely.We finally focus on the applicability of min/max selectors as the logic elements to switch between active constraint regions, which were left undetermined in Fig. 3 as ''select'' blocks.These selectors are applied on the controller outputs    and  0  associated with the controlled variables   and  0  , respectively, resulting in the process input (MV)   to be applied to the system.This methodology was Fig. 3. Decentralized control structure for optimal operation according to Theorem 2. The ''select'' blocks are usually max or min selectors (see Theorem 3).adopted in Krishnamoorthy and Skogestad [16] for optimal operation in the scalar case, ..with a single MV, where it was concluded that a constraint with a positive gain (  > 0) requires a min selector, whereas a negative gain (  < 0) requires a max selector.In the next theorem, we present similar results for the multivariable case.(𝑑), with   > 0 and arbitrary  , (), and that   is constant.Consider the control structure in Fig. 3 and Theorem 2, and assume that every possible subsystem is stable.

Theorem 3 (Decentralized Control. Applicability of Min/Max Selectors). In addition to Assumption 1, assume that the Hessian of the cost function with respect to the inputs is constant and positive definite, that is,
Let  0  denote the value of   that controls  0  =    ∇   (, ) = 0, and let    denote the value of   that controls   (, ) = 0.For a given active set , the associated nullspace of the active constraint gain matrix . Define the scaled projection matrix   and the transformed constraint gain matrix     as: The optimal input is given by  *  = min ( 0  ,    ) if the th diagonal element of the transformed gain matrix is positive ((    )  > 0) for any active set  that does not include .Conversely, the optimal input is given by Proof.See Appendix A □ )  changes sign for different active sets, a single type of selector would not account for all theoretical regions.The singleinput case of Theorem 3 can be easily verified by writing   −  0 = − 1       .As   > 0 for a convex optimization problem,   > 0 leads to  * = min ( 0 ,   ), and   < 0 leads to  * = max ( 0 ,   ), which is equivalent to the result in Krishnamoorthy and Skogestad [16].
Can some of the assumptions in Theorem 3 be removed?According to Theorem 3, the use of max-(or min-) selectors in Fig. 3 assumes that (    )  remains positive (or negative) for any active set  that does not include .This is to rule out cases where the steady-state gain for control of the constraint   changes sign, as this would lead to instability with integral action in the controller.In other words, this is to rule out interacting processes where  *  = min ( 0  ,    ) for a given active set , and  *  = max ( 0  ,    ) for another.However, it is not clear whether this is a restriction in practice.Thus, it is possible that the assumption about no sign change for the diagonal elements (    )  is not needed.This is left as an open research issue.
Cascade implementation.It is anyway possible to avoid this restriction by using the cascade switching implementation in Fig. 4.That is, for this implementation the simple selector logic is always optimal without the assumption about the sign of (    )  in Theorem 3. In the cascade implementation in Fig. 4, the constraints are always controlled in the lower layer, and the optimal constraint setpoint    will either be the value  ,0  that controls  0  = 0 or the constraint's limit value itself, such that    = min ( ,0  , 0) leads to optimal operation.This result can also be obtained by rewriting Theorem 3 in terms of , where it can be verified that the condition ( ,   )  > 0 is always satisfied.This result is presented in Appendix B.
The idea of using cascade control for self-optimizing control and constraint satisfaction has been previously proposed in Cao [7].There, the cost gradient is controlled in the unconstrained case, while the lower layer keeps the system feasible by saturating the setpoint from the upper layer.This approach ensures feasibility and self-optimizing behavior at the unconstrained region, but optimality at all active constraint regions is only ensured by carefully selecting the CVs at the upper layer, which is the main idea of the present work.In addition, even though the cascade structure will always operate optimally at steady state, it requires that the outer controllers  0   are sufficiently slower than the inner controllers    , and therefore the generic structure in Fig. 3 offers more flexibility in terms of loop tuning and implementation of further control overrides.The simulations presented in this paper are for the implementation in Fig. 3.

Case study 1 -Toy example
In this section, to illustrate the implementation of the proposed control structure, we consider a linear process with a quadratic cost function and 2 linear constraints.The process has 2 dynamic states , 3 inputs (MVs) , and 2 disturbances .The linear state-space model is: with  1 = 1 and  2 = 2.It is assumed that both states are measured, that is,  =  +  with  =  and  = 0.The steady-state optimization problem in terms of the states is: At steady state, the states can be eliminated to give the following static optimization problem: For given disturbances , we can solve the problem in Eq. ( 15) to find the optimal steady-state inputs  * and the active set . From this, we can graphically represent the active constraint regions as a function of the two disturbances as shown in Fig. 5.Note that this is done for visualization purposes only and is not a part of the proposed method.In fact, for the proposed method we do not need to know what the disturbances are; what is needed is measured or estimated values for the constraints  and the unconstrained cost gradient ∇  Ĵ .We see in Fig. 5 that all 2   = 2 2 = 4 combinations of constraints are possible.Each region has a specific set of CVs for optimal operation, namely the active constraints and the corresponding reduced gradients, as given in Theorem 2. We have   = 3 and   = 2, so with the proposed method, we need to design   +   = 5 SISO controllers with   = 2 selectors to obtain optimal steady-state operation.Since   >   , we always have   −   = 1 unconstrained degree of freedom corresponding to the controlled variable  0 =  0  ∇  Ĵ .From the nullspace of the full   matrix, we find that this direction is given by In addition, there are two unconstrained directions related to the two constraints.We have that  0 1 =   1 ∇  Ĵ should be controlled when  1 is not active, and  0 2 =   2 ∇  Ĵ should be controlled when  2 is not active.These directions are: For designing a decentralized control structure, a pairing between the constraints and the MVs must be performed.From the steadystate gain matrix   we see that  3 should not be used to control  1 (because of zero gain).Otherwise, there are no clear restrictions, and  1 is arbitrarily paired to  1 , and  2 is paired to  2 .We must require that the corresponding unconstrained optimal CVs are paired accordingly, meaning that  0  1 is paired to  1 ,  0 2 is paired to  2 , and  0 is paired to  3 .
For selector design, Table 1 shows the transformed constraint gains calculated using Eq. ( 12) for all active constraint sets, and we verify that the gains are always positive for both constraints.This means that selectors are possible for both control loops and that both selectors should be ''min''-selectors.The resulting control structure is shown in Fig. 6.
The cost gradient is estimated through a relinearization of the dynamic model at the current estimated state to obtain the following linear model: where, by setting ẋ = 0, the estimated steady-state cost gradient becomes [13] ∇  Ĵ = −   −1  +   (17) Fig. 6.Decentralized control structure for case study 1.
For state estimation, the model is augmented to include the disturbances as integrating states, according to: To estimate the states, a continuous-time Kalman filter is implemented with this augmented model, and the estimated state x and current input  are used to evaluate the matrices in Eq. ( 16) at all times, leading to the estimated cost gradient ∇  Ĵ in (17) .The matrices , ,   , , and  are as defined in Eq. ( 13), and   and   are calculated from Eq. ( 14) to give: We emphasize that analytical expressions for these derivatives are available due to the simplicity of this case study, and we encourage the use of automatic differentiation tools to obtain these matrices in more realistic case studies.
The constraint controllers   1 and   2 were designed according to the SIMC rules [18] with the choice  ,1 = 0.1 s and  ,2 = 0.01 s.In terms of gradient control, we assume that the effect of the inputs on the estimated ∇  Ĵ is that of a pure gain process, neglecting any dynamics associated with the gradient estimation, and therefore the gradient controllers  0  1 ,  0 2 , and  0 become pure integral controllers.These were tuned according to the SIMC rules with   = 0.5 s.All controllers linked to selectors are implemented with anti-windup action  based on the back-calculation strategy [14,19], with a tracking time of   = 0.01 s.The resulting controller tunings are summarized in Table 2.
The closed-loop simulations are shown in Fig. 7. To validate the optimality of the control structure, the disturbances were changed stepwise every 15 s (see lower left plots) to make the system operate in all four active constraint regions (see lower right plot).It can be seen that constraint changes are effectively handled, giving up the corresponding gradient projection when a constraint becomes active, and that operation is driven to the optimal steady state for all disturbances.

Case study 2 -Williams-Otto reactor
The control structure proposed in Section 2 depends on using projection matrices.These are constant only when the constraints are linear in the MVs.We now consider a nonlinear case study where this assumption is not satisfied and one may expect economic losses in some regions.The case study is based on the process described by Williams and Otto [20] and studied in [16], see Fig. 8.It consists of a continuously stirred reactor tank with perfect level control, in which A and B are mixed, generating the main product P, the less interesting product E and the undesired byproduct G.The reactions and reaction rates are given by: The component mass balances for the six components give the following set of ODEs:  ]  , where   is the relative change in the price of the main product, as defined in Eq. (20).
The active constraint regions as a function of the two disturbances are shown in Fig. 9.In contrast to the previous case study, the lines delimiting each region are not straight.This alone should not affect the optimality of the proposed framework, as the optimality only requires that the constraints are linear in the MVs.However, since the latter does not hold for the case study, the use of constant projection matrices will lead to some economic loss.
We have   =   = 2 so Assumption 1 is satisfied.With the proposed method we need to design   +   = 4 SISO controllers with   = 2   Since   =   , the system has no completely unconstrained degrees of freedom, so there are no variables  0 that are always controlled.The gradient projections  1 and  2 for the two potentially unconstrained degrees of freedom become: For MV-CV pairing, we choose  1 =   for controlling  1 and  0 1 =   1 ∇  Ĵ , and  2 =   controlling  2 and  0 2 =   2 ∇  Ĵ .This pairing choice was made based on the steady-state RGA for constraint control, which gives  = 0.638 for the chosen pairing.For designing the selectors according to Theorem 3, a local analysis of the transformed constraint gains given in Eq. ( 12) was made at the nominal point and is summarized in Table 5.The projected gains are negative for both constraints regardless of the active set, which means that both selectors should be ''max'' selectors.The resulting control structure is presented in Fig. 10.
To tune the controllers, we obtained the following transfer functions from the MVs to the constraints (with time in hours): ,  22 (s) = −0.002410.072 s + 1 Based on this, the PI controllers for the constraints were tuned using the SIMC rules [18] with  ,1 = 0.005 h and  ,2 = 0.01 h.Similar to the previous case study, the method for gradient estimation is again considered to be instantaneous with respect to the inputs, meaning that the gradient controllers become integral controllers, tuned using the SIMC rules with   = 0.05 h.All four controllers were implemented with anti-windup action based on back-calculation with a tracking time of   = 0.01 h.The controller tunings are summarized in Table 6.
Closed-loop dynamic simulations are presented in Fig. 11.The disturbances were changed so that all four active constraint regions were explored.Since we consider that the cost gradient ∇  Ĵ is an available measurement (from an estimator with a perfect model), operation in the fully unconstrained region (from  = 21 h to  = 27 h) is optimal at steady state, which can be seen by the input values converging to the exact steady-state optimal value.Since we assume that the constraints are directly measured (which is a mild assumption), the same logic applies to the fully constrained region from  = 9 h to  = 15 h.In addition, operation is optimal at the nominal point by design.In the two remaining partly constrained regions, the system does not converge exactly to the steady-state optimum, but the constraints are always satisfied (except for short dynamic transients, which may be avoided by introducing a back-off for the constraints).It is interesting to note that for the third set of disturbances ( = [1.0,−0.2]  from  = 6 h to 9 h), the second constraint ( 2 = 0) is not controlled, even though it should be optimally controlled together with  1 = 0. Instead, the selector logic results in the control of  0 2 = 0, which can be done without violation of  2 , that is, constraint  2 is ''over-satisfied''.The reason for this non-optimal operation is that the selected value for projection matrix  2 is not optimal in this operating region.
The steady-state economic loss is better visualized as a function of the disturbances in Fig. 12.The highest losses are observed around where we ideally should switch between the partly constrained and the fully constrained regions.The optimal switch between these regions (black lines) does not coincide with the actual switch obtained with the selectors (blue lines).Economic loss is observed before the optimal switch due to the inaccuracy of the projection matrices.For the same reason, and because this further leads to suboptimal performance of the selectors, economic loss is also seen between the optimal and actual switch of CVs.However, the optimal switch between the fully unconstrained and the partly constrained regions coincides with the actual switch between the corresponding CVs.This happens because, before the switch, the full cost gradient ∇  Ĵ is controlled to zero, leading to zero economic loss, and the constraint becomes active immediately at the switch.Therefore, at this switch, the economic loss is zero, and it continuously grows as the system moves further into the partly constrained region.

Steady-state cost gradient estimation
The results in this paper assume the availability of the steadystate cost gradient ∇   during operation.This can be fulfilled through model-based estimation, model-free estimation, or a combination of both methods [3].Model-free methods usually depend on the perturbation of the inputs, and when the constraints are being controlled the perturbation can be done in their setpoints instead.In the presented case studies, we used a model-based approach, where a Kalman Filter was used to estimate the current dynamic state  and disturbance  with an augmented model ( 18), and then setting ẋ = 0 in the linearized model ( 16) leads to the gradient estimate (17).
Because the cost gradient ∇   is, by definition, a steady-state variable, it is not well-defined during a dynamic transition, and any gradient estimator must make some steady-state assumption or prediction.The necessity of estimating the cost gradient is related to ensuring exact optimality.In practice, one would wish to use an approximation of the cost gradient that is more easily implementable, even if that means accepting some economic loss.In that sense, data-driven approaches for this estimation would be appealing, as well as self-optimizing control methods that provide an approximation for the cost gradient through a static combination of measurements [6].
However, a simpler approach is to use a static estimation of ∇   directly based on the measurements .In another paper [21], we prove the optimality of a simple linear steady-state gradient estimate of the form where  are the measurements, and the constant vector   and the constant matrix  are obtained using the ''exact local method'' of selfoptimizing control.In addition, a correction of   from a more accurate gradient estimator may be applied on a slower time scale, for example, using a model-based approach like RTO or a data-based perturbation method like extremum-seeking control.

Handling of constraints
Constraints in process systems are usually measured or estimated, and our approach is optimal for such cases, as control loops are implemented to handle these constraints.In MPC applications, where the problem is formulated as a dynamic trajectory optimization, it is common that process constraints are posed as constraints on the dynamic states, but this is not how our approach handles this issue.Rather, our method is focused on process constraints that may become active at steady state and influence process economics.

Updating of projection matrices
In the simulations presented in this paper, we assumed that a linearization of the constraints at a nominal operating point would be sufficiently accurate for capturing the transitions between active constraint regions.This simplification was primarily made to ensure a control structure that can be easily implemented, and it led to acceptable results even for a nonlinear case study (see Fig. 12).However, it is possible to enhance economic performance by updating the projection matrices  and  0 during operation.To accomplish this, an accurate estimator for the complete constraint gradient matrix   is required, and typically such an estimator is only available at a time scale similar to that of RTO.However, our primary objective is to achieve acceptable economic losses in fast time scales, and this could be accomplished by using constant projection matrices.

Controller tuning
Even though the proposed control structure only has   +   SISO controllers to be designed, the tuning of these controllers may prove to be challenging.This is because these controllers must work in many different regions (up to 2   theoretical regions), and the interaction between loops will change depending on which controllers are active.The pairing between MVs and CVs should consider this, and the tuning for the loops should be robust in the sense that acceptable performance is attained for every operation mode.This issue was not noticed in the case studies in this work, but it is easy to see that it may arise in practice.

Limitations for systems with many constraints
In this paper, we consider a class of problems with   ≥   , so it is possible to devise a simple, decentralized control structure.There is a particular case of systems with more constraints than inputs that can fit into the framework proposed in this work.That would be the case where the constraints can be arranged into   groups, where each group is comprised of constraints that have parallel gain vectors respect to the inputs, .. the constraints   and   would belong to the same group if ∇    = ∇    for some nonzero .In practice, this would represent a process variable with lower and upper bounds, or constraints of similar nature caused by different factors, .. a maximum processing rate due to upstream or downstream conditions.Each of these groups has a unique characteristic direction in terms of the rows of   , which can be used to calculate the gradient projections with the methodology described in Theorem 2. Each of these groups should then be organized internally following the single-variable methodology described by Krishnamoorthy and Skogestad [16].As the methodology devised in this paper mitigates the correlation between each group of constraints, the gradient projections that serve as unconstrained CVs remain constant with respect to changes in the remaining loops, and therefore no additional logic is required in the implementation of max/min selectors for constraint handling.
The main case not covered by the present methodology is when there are   >   independent constraints that may become active, expressly violating Assumption 1.In this case, considering that some pairing between MVs and constraints is done, the first problem that arises is the possibility of constraints paired with the same MV becoming active at the same time, requiring that one or several constraints become controlled by other MVs.A heat exchanger case with   = 3 >   = 2 was studied in Bernardino et al. [22], where it was shown that a region-based approach similar to the one studied in the present paper fails to achieve optimality for some disturbance scenarios, whereas the primal-dual approach always reaches the optimal steady state.To achieve optimality with a region-based approach, an adaptive pairing strategy may be used, as described for this case study in Bernardino et al. [23].This gives optimal operation for all disturbances, but the adaptive pairing becomes quite complicated (see Figure 2 in [23]).
In this sense, general strategies for switching pairings require more complex logic, and currently, there is no systematic arranging of classical control logic blocks that can account for that.On top of that, even if conflicting constraints are not an issue for the considered operating window, ..constraints paired to the same MV do not become active at the same time for the considered disturbances, the design of controllers for the unconstrained degrees of freedom becomes more complicated.As the constraints are assumed to be independent, the gradient projections optimally controlled to zero will be different when each of them is active.This entails that the remaining control loops have to change depending on which controller related to this MV is active.Therefore, proposing decentralized control structures for the optimal operation of systems with more constraints than MVs inevitably leads to complex and interacting control loops, and centralized strategies such as the primal-dual feedback optimizing control presented in Krishnamoorthy [24] or MPC become more appealing.
For the same reason, the proposed framework has limitations in optimal dealing with input saturation.In real systems, every MV has physical bounds   ≤  ≤   in addition to the process constraints .Therefore, every physical system in a way has more constraints than MVs, and one must identify the constraints that are more likely to become active if one follows Theorem 2 for designing a control structure.The choice of not pairing an MV that may saturate with an important CV, in this case, an economic constraint, agrees with the rule of thumb ''pair an MV that may saturate with a CV that may be given up'' [25], as the gradient projections paired to that MV should by design be given up in case of MV saturation.
The proposed framework attains optimal operation in a wide operating range, by enforcing optimality conditions at steady state for all possible active set combinations for a maximum of   independent constraints.It should also be emphasized that less frequent constraints can still be dealt with in the current framework by the implementation of more selectors, even if   >   , bearing in mind that steady-state optimal operation will not be guaranteed when those become active due to changes in the unconstrained CVs.However, violation of such infrequent constraints would be prevented, which is the main goal of such additional control loops.

Stability and optimal convergence of selectors
The control strategy proposed in this work relies on switching blocks to perform optimal operation in different operating regions.Analyzing the stability of switched systems is more complex, as the stability of a switched system may not necessarily match that of its corresponding continuous subsystems [26].In Theorem 3, we assume that each subsystem within the switching system is stable, which is a condition already present and well described when using decentralized control in multivariable systems.By ensuring that every subsystem is stable, the overall stability of the switching system can be guaranteed.This can be achieved by enforcing a sufficiently large average dwell time [27], which is a practical and easily implementable solution.
The conditions for implementing min/max selectors to detect switches in active constraints optimally are outlined in Theorem 3.This theorem is based on a local analysis of the optimization problem and is rigorously applicable to problems with a constant positive definite Hessian   and constant constraints gain   .A relevant case in practice is that of linear economic objectives, for which   is positive semidefinite, but this case is always solved by active constraint control, as there are no unconstrained degrees of freedom to be determined.While the presented proof does not address generic nonlinear optimization problems, it provides a useful local test that can eliminate certain impossible configurations resulting from the chosen MV-CV pairings or the formulation of the optimal operation problem itself.If the conditions specified in Theorem 3 are not satisfied, we recommend utilizing the cascade framework presented in Fig. 4.
The condition derived in Theorem 3 for applicability of selectors would only be violated by highly interacting systems, where the sign of the transformed constraint gain (    ) would change depending on the active loops.This condition is conjectured to be associated with the decentralized integral controllability (DIC) of each potential subsystem [28].In our study, we could not find an example of a linear system with a convex objective function and without DIC that does not satisfy the conditions stated in the selector theorem.This further suggests a connection between these concepts and that the DIC conditions possibly imply the applicability of selectors.The link between the conditions of Theorem 3 and controllability aspects remains an open challenge that requires further investigation.

Other switching approaches
In this work, we have proposed the use of selectors in the controller outputs for detecting switches in active constraints.However, other strategies for adaptively controlling constraints in the context of optimal operation have been proposed.Manum and Skogestad [10] studied the problem of active constraint switching in self-optimizing control by tracking the self-optimizing CVs in neighboring regions, where the switching happens when there is a change of sign in the monitored variable.In the notation herein presented, this would be equivalent to the following switching logic: The problem with implementing such logic lies in the resulting dynamics of the control system.As this logic implies that the reference variable is perfectly controlled for accurate detection, the logic should operate in a slower time scale than that of the closed-loop system, which would in turn result in undesired behavior, especially constraint violation.Operating the switching logic in fast time scales could in turn lead to the appearance of limit cycles, due to self-sustained switching between control loops.
We have also presented the cascade control structure in Fig. 4 as an alternative switching strategy.A similar idea has been proposed by Cao [7] to promote self-optimizing operation at the nominal region while sub-optimally coping with constraint satisfaction.There are however some disadvantages to this approach related to the limitations that the cascade structure imposes.If constraint control is slow, controlling the corresponding gradient projection becomes unnecessarily slow.Moreover, even though constraint control may help with decoupling the system, it may also cause the opposite problem, and the interaction between loops may impose limitations on the performance of the upper layer.Therefore, the use of a cascade framework for optimal operation may be beneficial, but the improvement that it may bring must be assessed for each particular case study.
Recently, the work of Ye et al. [29] has tackled the problem of changing active constraints by embedding the switching constraints into the CV design, generating a single nonlinear CV.Because the resulting CV design problem was deemed intractable in most cases, a neural network was used to approximate the behavior of this theoretical CV.It is interesting to note that the switching behavior still happens in the designed CV, with the exact ideal CV being in general non-smooth.This is expected because of the nature of the problem, and although neural networks can approximate these variables, the interpretability of the resulting CV is lost, and constraint control must be explicitly performed elsewhere.In the present work, we deal with the switching explicitly, controlling the constraints directly when it is optimal.

Conclusion
We propose a simple framework for decentralized optimizing control with changing active constraints.The starting point is that at steady state, optimal economic operation in a given active constraint region  is achieved by keeping the controlled variables  = [   ;    ∇   ( * , ) ] at constant setpoints   = 0 (Theorem 1 [6]).Here   denotes the set of active steady-state constraints, and    ∇   ( * , ) is the reduced steady-state cost gradient for the remaining unconstrained degrees of freedom.
There are some degrees of freedom in the choice of the directions in the unconstrained nullspace   and to implement constraint switching in a simple manner, these should be chosen in accordance with the constraint directions.The main contribution of this paper is to prove in Theorem 2 for the case with   ≥   , that we should control the unconstrained variables  0 =  0  ∇   (which are not affected by the constraints), and in addition, depending on whether the constraint   is active or not, either control the constraint    =   or the associated unconstrained variable  0  =    ∇   , where  0 and   are calculated according to ( 7)-( 9).This can be implemented with the simple control structure in Fig. 2. Furthermore, Theorem 3 shows that the switching can be performed with min/max selectors, which leads to the simple control structure in Fig. 3. Here, no centralized supervisor is needed to determine the active constraints, as the switching logic uses local feedback controllers.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
L.F. Bernardino and S. Skogestad Substituting this equation into the second optimality condition gives the following expression for   * : The variables can be partitioned as follows: ) is equivalent to a scaled projection to the nullspace of    ,  () [17], according to the identity   =  () (   ()    () , and is therefore positive semidefinite.The effect of inclusion of an arbitrary constraint   is therefore given by: We can see that the th component of the vector   (   )  dictates the steady-state behavior of the th MV when the th constraint becomes active and the system is operating at the active set .Following the notation introduced in the statement of Theorem 3 and in Fig. 3, we have    = (  * )  and  0  = (  )  , since these are the MV values such that   = 0 and  0  =    ∇   = 0, respectively.This means that, for (  (   )  )  > 0,    −  0  < 0 when   > 0 and consequently    should be selected, and    −  0  > 0 when   < 0 and consequently  0  should be selected, meaning that we must use  *  = min ( 0  ,    ) for guaranteeing optimality.Similar analysis can be performed for (  (   )  )  < 0, leading to  *  = max ( 0  ,    ).Since this analysis is performed for arbitrary  and  ∌ , guaranteeing that (  (   )  )  = (    )  has the same sign for any possible  is sufficient for guaranteeing the theorem statement, which completes the proof.□

Appendix B. Optimality of min selectors in cascade structure
In this section, we restate Theorem 3 for the cascade case illustrated in Fig. 4, proving its optimality.In this control structure, the manipulated variables as seen from the higher layer can be represented at steady-state as: Note that we must require that  , the Jacobian for the change of variables, is full rank, such that the optimality conditions for Eq.(B.1) and Eq. ( 1) are equivalent.This is a mild assumption related to the steady-state controllability of the constraints in the lower layer with the chosen pairing, and it results in a transformed Hessian   =      that is also positive definite.Also, for the transformed problem, we have the gain matrix  , from the transformed inputs to the constraints: This allows us to write the optimal CVs in terms of the projection matrices   and  0 for the transformed problem as: The procedure for obtaining the difference in the optimal solution for neighboring regions presented in Appendix A is also valid for the transformed problem, and we must therefore analyze the sign of the diagonal of the matrix product  ,   .Recall that   =  () (   ()    () ) −1   () , where the matrix   ()    () is positive definite, being here a principal submatrix of   which selects the rows and columns with indexes not in .Therefore,   becomes a positive semidefinite   ×   matrix, where the diagonal elements    are zero for  ∈  and positive for  ∉ .
Finally, we can see that ( ,   )  > 0 for  ∉ , since the first   elements of the diagonal of  ,   are the same as those of   .It follows that, according to Eq. (A.1) and the results here obtained, the optimal solution is given by  *  = min ( 0  , 0) for  = 1, … ,   , which completes the proof.

Fig. 4 .
Fig. 4. Decentralized control structure for optimal operation, using an alternative cascade implementation.

Fig. 5 .
Fig. 5. Active constraint regions for case study 1 as a function of disturbances.
L.F.Bernardino and S. Skogestad

Fig. 9 .
Fig. 9. Active constraint regions for case study 2 as a function of disturbances.

11 .
Closed-loop simulation results for case study 2.

12 .
Steady-state closed-loop economic loss for case study 2.

Table 1
Diagonal elements of

Table 2
PI controller tunings for case study 1.Note that   =   ∕  is the integral gain.The first four controllers have anti-windup with tracking time   = 0.01 s.

Table 3 .
The economic cost  includes the cost of reactants   and   and the selling price of products   and   , and the operational constraints are related to maximum allowed values for   and   .The steady-state optimization problem becomes min     +     − (  +   ) [   (1 +   )  +     ] =

Table 4
Nominal operating point for case study 2.
steady-state model of the constraints.In the following simulations, we use the linearization performed at the nominal operating point presented in Table4, leading to fixed CVs for operation in all regions.

Table 5
Diagonal of

Table 6
PI controller tunings for case study 2.
is being controlled to zero, a change of sign in   means that the th constraint became active, as this sign change corresponds to constraint violation; • Conversely, if   is being controlled to zero, a change of sign in  0  =    ∇   means that the th constraint became inactive, as this sign change corresponds to a change in the objective function slope.