Decentralized PID Controller Tuning Based on Nonlinear Optimization to Minimize the Disturbance Effects in Coupled Loops

Decentralized proportional-integral-derivative (PID) control systems are widely used for multiple-input multiple-output (MIMO) control problems. However, decentralized controllers cannot suppress the plant interactions in multivariable systems, which are addressed in the controller tuning phase. In this paper, a decentralized PID tuning method is proposed in order to minimize the undesirable effects of the coupling between the inputs and outputs of the closed-loop system. For this purpose, the PID parameter tuning method solves a nonlinear optimization problem. This optimization problem is formulated with the criteria of the performance, robustness and multivariable stability of the closed-loop system. A single design parameter is required to specify the trade-off between performance and robustness. Simulation studies are conducted to demonstrate the effectiveness of the proposed method. The performance is compared to that of alternative tuning techniques from the literature. Results show that the proposed approach is a feasible candidate for industrial application, as it is simple to implement and capable of addressing robustness and stability concerns of plant operators.


I. INTRODUCTION
Decentralized proportional-integral-derivative (PID) controllers are frequently used in the regulatory control layer of industrial process plants. Single-input single-output (SISO) control is the dominating control structure in the industry for many reasons: it maintains the simplicity of the control system; it is easier to maintain; there are fewer tuning parameters compared to full multivariable controllers; it provides flexibility; and it can easily be made fault-tolerant [1]. Moreover, even when multiple-input multiple-output (MIMO) control strategies are used, such as model predictive control (MPC), they typically operate in a supervisory mode with decentralized PID controllers at the lower level [2].
Regardless of their practical benefits, single-loop controllers cannot suppress the interactions between variables in MIMO processes; thus, every system input affects every system output to some extent, unless static and dynamic The associate editor coordinating the review of this manuscript and approving it for publication was Qi Zhou. relationships between a given input-output pair are nonexistent. Since the SISO controllers interact with each other, tuning each loop independently is not recommended. Applying SISO tuning methods for a MIMO system often leads to poor performance and stability [3].
In this paper, a nonlinear constrained optimization problem is proposed to compute decentralized PID parameters, considering the trade-off between performance and robustness. The proposed problem consists of minimizing the disturbance effects in the output of each loop due to changes in the operating point of coupled loops. For this purpose, the chosen cost function is the sum of the integral of the absolute error (IAE) from each output during a disturbance due to step set-point changes. The robustness to modeling errors for each loop is considered based on a constraint on the peak of the sensibility function. Moreover, the overall stability is guaranteed by a constraint based on the biggest log-modulus criterion. Although we cannot guarantee that it finds the globally optimal controller parameter values, through many examples, the proposed method always meets the robustness VOLUME 9, 2021 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ requirements specified by the designer. To the best of our knowledge, no studies in the literature combine these three criteria into a single optimization problem to tune decentralized PID controllers. Some benefits of the proposed method, compared to other approaches from the literature, are as follows: (i) this method is not restricted to diagonally dominant processes, (ii) model reduction is not necessary, (iii) high-dimensional processes are not an issue, (iv) all controller parameters are computed at once with guaranteed overall stability, (v) no detuning factor is necessary, and (vi) a robustness and performance trade-off is specified for each loop by only one tuning parameter, which varies based on the operating settings of the process.
The remainder of this paper is structured as follows. In Section II, a literature review of the main contributions to decentralized PID tuning using optimization is discussed. Section III outlines all assumptions regarding the multivariable system and the control system structure. Section IV introduces the performance and robustness criteria used together in the constrained optimization problem stated in Section V. Simulation results are presented in Section VI, followed by the conclusions in Section VII.

II. LITERATURE REVIEW
Several researchers have paid attention to the tuning of decentralized PI/PID controllers over the years. The biggest logmodulus tuning (BLT) method [4] is a well-known method that considers the interactions among the loops. The BLT method designs controllers by applying the Ziegler-Nichols method to the diagonal transfer functions and then detunes them by introducing a single detuning parameter to meet the stability criterion of the biggest log modulus. This method provides reasonable preliminary controller settings with guaranteed closed-loop stability. Chen and Seborg modified the Ziegler-Nichols tuning rules based on the ultimate point for a Gershgorin band in [5], but their method still uses a detuning parameter, similar to the BLT method. [6] developed an iterative method for decentralized PID controller tuning based on the effective open-loop transmission from input u i to output y i and phase/gain margin specifications. The iterative procedure solves only one optimization problem for each controller and is applied to all loops, and its solution is obtained by linear programming. In the same vein, [7] proposes an iterative design methodology of multiloop PID controllers. The proposed method uses a frequency response matrix representation of the system to avoid process approximations. Furthermore, different frequency domain robustness margins are used as specifications. The EOP approach was also used in [8] for high dimensional MIMO systems but was restricted to PI controllers. Alternatively, to reduce the effect of process interactions in a decentralized control system, [9] proposed the use of the double-loop multiple scale control scheme (DL-MSC).
Several studies from the literature have cast the PID tuning problem as an optimization problem in order to improve the performance of the closed-loop system by means of decentralized PID tuning. In [10], the concept of Gershgorin bands was used, and a nonlinear optimization problem was formulated with constraints on stability margins, but the method was limited to weakly coupled processes. In [11], a nonlinear optimization problem was formulated, but the constraints were related to performance indices such as the maximum overshoot and maximum controller limit deviation for each loop; however, no overall stability was considered for the MIMO closed-loop system. In [12], two linear programming approaches were proposed to compute the PID parameters, where the loop interactions were taken into account by the Gershgorin bands and the effective open-loop process concept, but each loop had to be tuned independently.
Centralized PID approaches are also present in the literature, and in contrast to the previously discussed decentralized approaches, they are often realized by means of a decoupler or some other mathematical tool that tries to mitigate variable coupling effects. However, a decoupler makes the structure of the control system more complex and highly dependent on the plant models due to the need to know the process model to design the decoupler. [13] studied centralized PID control of TITO (Two-Input Two-Output) systems. The method proposed by the authors designs PI controllers for the main diagonal based on [14] and decouplers for the off-diagonal dynamics based on model approximations. The combination of both effects can successfully control a coupled tank system, delivering performance and robustness in servo and regulatory operation. [15] extended the BLT method to centralized control of MIMO systems. In [16], the resulting controller is based on the internal model control (IMC) structure. The process is identified as a FOPDT transfer function and can account for disturbances and model uncertainty. In contrast, [17] presents a novel centralized controller to MIMO industrial processes with heavy interactions and significant time-delays. The main drawback of implementing a centralized PID approach is that controller stability is only guaranteed if all of the loops are permanently closed. In industrial practice, this is a weak assumption since keeping some loops in manual mode or out of service is common due to maintenance procedures, while other loops may be in regular operation. In this paper, only the design of decentralized multivariable control systems is considered.
In the vein of multiobjective optimization, [1] utilized an evolutionary algorithm (bat algorithm) to design a fractal PI/PID controller for TITO. In [18], the design and tuning of a PI controller for a nonlinear TITO process based on different optimization approaches, such as a genetic algorithm (GA) and particle swarm optimization (PSO), was investigated. The authors demonstrated that for the nonlinear case, the GA yields faster tracking performance and less overshoot than the PSO-based solution. This result might indicate that a GA is preferred for tuning purposes in the presence of high nonlinearities and model mismatch. In contrast, [19] presented an optimization method of tuning decentralized PI/PID controllers based on GA. Simulation results demonstrate that the decentralized PI control was compatible to the referenced method while the decentralized PID control was better than the referenced method.
In [20], the authors recast the MIMO PID tuning problem as the calculation of an optimal detuning factor by means of a differential evolution algorithm. This strategy considerably reduces the search space, making the method faster, especially for larger order systems. The authors illustrated the technique on a 4 × 4 system, highlighting its computational efficiency, which can lead to successful applications in complex and time-sensitive industrial applications. The fractal order PID approach works better for some processes according to the authors. [21] utilized a GA to design fractal PID controllers for MIMO systems; however, the resulting optimization problem became computationally expensive due to the high number of decision variables. Meanwhile, [22] utilized a multiobjective variant of the PSO algorithm to address the same problem. Results showed that obtaining the set of Pareto optimal solutions is difficult for such a problem. Also based on optimization problems, [23] proposed a method to design decentralized and centralized PID controllers. In this case, the objective is to minimize the low-frequency sensitivity of the closed-loop system using LMI. [24] developed a tuning technique for centralized and decentralized PI controllers based on goal attainment. The centralized approach dramatically increases the closed-loop performance in an open-loop unstable TITO system, but the drawback of the goal attainment method is that it requires a priori definition of the goals, and the final solution is heavily dependent on this definition.

Consider the MIMO system with n inputs (actuators) and n outputs (sensors) represented by the transfer function matrix
where G ij (s) represents the system dynamics between the input j and output i, which is assumed to be linear time invariant and strictly proper. The structure of the decentralized PID controller is given by The closed-loop system is shown in Figure 1. The vector r denotes the set-points, u is the vector of system inputs (manipulated variables), and y is the vector of process outputs (controlled variables). The control structure is assumed to have been correctly chosen in a previous design based on some measure of interaction, for instance, using the relative gain array (RGA) [25] or one of its variants. The PID controller considered in this paper is formulated as where k p j , k i j , and k d j are the proportional, integral and derivative gains, respectively, for each j controller, j = 1, . . . , n.
The constant T f j > 0 is the derivative action time constant and is assumed to be fixed. Although the PID formulation in (3) is considered in this paper, the proposed method is not restricted to this formulation. Therefore, different PID formulations can be directly used without requiring parameter conversions. This feature is particularly interesting for industrial practice since commercially available PID packages might take different parametrizations.
The problem can be stated as follows: Under the system structure shown in Figure 1, define the PID controller parameters in order to minimize the disturbance effects on the outputs due to set-point changes in coupled loops. For this, the trade-off between performance, robustness and multivariable stability are considered into a single optimization problem.

IV. EVALUATION CRITERIA
The proposed method is based on three criteria, each of which is related to the essential performance and robustness characteristics of the multivariable system. The aim is to reduce the disturbance effects on the outputs due to set-point changes in coupled loops. Each of the chosen criteria is presented in the following subsections.

A. THE PERFORMANCE CRITERIA
Due to interactions in the MIMO process, a change in the input of loop i affects not only its corresponding output but also the other outputs of the process. This change in loop i is considered as a load disturbance by the other loops. For multivariable systems, disturbance rejection due to coupling is one of the main concerns since the disturbances directly affect the performance of the system.
The integral of absolute error (IAE) index is commonly used to characterize the load disturbance response [26]. Here, the IAE index is used to quantify the disturbance effects on the loop outputs due to the set-point change in the coupled loop. In particular, the IAE is computed for the load disturbance effects on the i − th output caused by a unit step change VOLUME 9, 2021 in the j − th set-point, i = j. Thus, the suggested criterion can be defined as where e i (t) = r i (t) − y i (t). The variable t indicates the time.
To illustrate the disturbance effect due to the set-point change in the coupled loop, consider the Wood-Berry distillation column [27]. This process is given by The following PI controllers, proposed by [4], are used: As shown in Figure 2a, a unit step change is applied in the set-point of loop 1. Due to the transfer function G 21 , the change in loop 1 also affects the output of loop 2. However, the controller C 2 acts to reduce this disturbance. Similarly, the controller C 1 acts to reduce the effect of a disturbance in u 2 on y 1 due to a unit step change in the set-point of loop 2 (see Figure 2b).
Note that the IAE value, as defined in (4), is equal to the area under the curve of the load disturbance response, as illustrated in Figure 2. In this way, the IAE index can be used to measure the disturbance effects on the loop outputs due to the set-point change in the coupled loops. The more aggressive the PID controller tuning for changes in the set-point is, the higher the effects of the disturbances in the coupled loop. In contrast, smooth controller tuning produces a slow response to load disturbances. Therefore, for all closed loops, the trade-off between performance and robustness should always be considered.

B. THE ROBUSTNESS CRITERIA
Robustness is an important feature of the closed-loop system since industrial processes operate under a wide range of operating conditions. In such cases, the controller should be able to stabilize the system for all operating conditions [28].
To define the robustness criteria, consider the loop transfer function for the j − th loop given by The sensitivity function for the j − th loop is defined as The maximum sensitivity is then given by where ω is the frequency in [rad/s] and M S j is an upper bound specification.
The maximum peak of the sensitivity function provides a useful measure of the general robustness [29]. In this way, the robustness criterion is defined as the largest value of the sensitivity function and is applied as an optimization constraint. This criterion provides good performance in terms of robustness with respect to the uncertainties in both the process model and disturbance. As a geometric interpretation, the value of S j ∞ is the inverse of the shortest distance from the Nyquist curve of the loop transfer function to the critical point (−1, 0). Figure 3a illustrates the loop transfer function from the same process for different M S j values. As shown in the Nyquist plot, for M S 1 = 1.45, the loop transfer function L j (s) is further from the critical point (−1, 0) than that for M S 2 = 2.0. Note that the circle expands with decreasing S j ∞ , resulting in greater closed-loop robustness and smoother output response to changes in the set-point, as illustrated in Figure 3b. Thus, the robustness constraint is fulfilled if the open-loop Nyquist curve does not enter the M S j circle. The closer the curve is to the critical point, the more aggressive the control system response is. An overview of frequency domain robustness constraints is presented in [30].
From the restriction in (8), the robustness of the MIMO system is not guaranteed. The main reason for establishing the robustness criterion for the individual loops is that this criterion is used to design a faster or smoother closed-loop response depending on process requirements. Additionally, the robustness is guaranteed in the case that one or more loops are opened, which is a common practice in the industry, for instance, for loop maintenance.
The trade-off between performance and robustness varies depending on the control problem. Therefore, having a design parameter to change the properties of the closed-loop system is desirable. High values of M S j lead to decreased robustness. According to [26], typical values of M S j are in the range of 1.2 to 2.0. M S j values larger than 4 lead to poor performance as well as poor robustness. The selection of M S j should take into account the model uncertainty. If the model uncertainty is significant, then M S j should be small (conservative closed-loop response).
For the proposed design method, the upper bound M S j on the peak of the sensitivity function is defined as the single design parameter. Once the peak of M S j is defined, gain and phase margin specifications are unnecessary. For instance, M S j = 2 implies a gain margin of GM > 2 and a phase gain of PM > 30. From a practical point of view, the definition of only one value to indicate robustness for each loop might be more acceptable by process control engineers.

C. THE MULTIVARIABLE STABILITY CRITERIA
Guaranteed stability of the individual loops is not a sufficient condition to ensure the stability of the multivariable closedloop system. Interactions between the loops may produce undesirable control actions, which may lead to overall system instability [29].
To ensure system stability when all single loops are closed, consider the characteristic equation of the multivariable closed-loop system: To analyze the encirclements of the closed-loop process with respect to the critical point (−1, 0), one must subtract 1 from (9) and obtain the scalar function The function (10) can be plotted in the complex plane as a function of frequency. The closer W (jω) is to the critical point (−1, 0), the closer the multivariable system is to closed-loop instability. Therefore, is equivalent to the closed-loop servo transfer function for SISO loops. To achieve stability for the overall system, [4] proposed the biggest log modulus, which is based on a multivariable stability criterion. This stability criterion is defined as follows: .
The biggest log modulus is a measure of how far the system is from closed-loop instability. In [4], it was suggested that for equal to 2n (where n is the number of loops to be tuned), the system provides reasonable responses for set-point changes and load disturbances. In this paper, the criterion is applied as an optimization constraint, and its upper bound is 2n.

V. THE OPTIMIZATION PROBLEM
Taking all three criteria described in Section IV, the constrained optimization problem is formulated as follows: which can be used to design all of the controllers C j (s), j = 1, . . . , n. Note that the only design parameter in the proposed technique is M S j , the upper bound on the peak of the sensitivity function for each j − th loop. The solution of the optimization problem in (13) is obtained by using the active-set method (see [31] for details), VOLUME 9, 2021 which is implemented using the fmincon routine from the MATLAB R optimization toolbox.
Due to the chosen cost function, the solution of the proposed optimization problem is not trivial. Its degree of difficulty depends on the number of control loops involved, which defines the number of decision variables of the problem. Additionally, high PID gain values may yield the fastest closed-loop tracking responses, but they may also incur large disturbances in the coupled loops. On the other hand, low PID gain values may lead to a slow response to load disturbances. Thus, the optimal solution is the one that better addresses the trade-off between performance and robustness considering the process as a whole.
The proper definition of the initial controller parameters is essential for a fast convergence of the proposed optimization problem. This initial controller ensures that the closed-loop system is in the stability region. For this, we recommend the use of the method formulated by [32]. However, other methods of tuning PID controllers can be used.

VI. SIMULATION RESULTS
In this section, three examples are considered to demonstrate the performance of the proposed method compared with those of other well-known methods. To ensure a fair comparison, the performance and robustness of the control system are measured using the evaluation criteria described in Section IV.
The examples share some common features. A grid of 1000 logarithmically spaced frequencies between ω min = 10 −2 [rad/s] and ω max = 10 2 [rad/s] is used in all examples. The initial guesses of the controller parameters are defined using the method proposed in [32]. For this, only the SISO model of the main diagonal of the multivariable process is considered.
In all cases, the IAE due to a unit step disturbance in the set-point is computed considering the entire simulation time. Additionally, the performance of the controllers is assessed using the cost criterion given by Another performance index used in all cases to evaluate the designed controllers is the sum of the absolute input increments (SII) of input i over the simulation horizon T , calculated as The greater the SII index, the more abrupt the changes in the variable manipulated by the controller are. An aggressive change of the manipulated variable may impair the lifetime of the actuator, and the actuator cannot reach such rapid variations in some cases.

A. EXAMPLE 1
Consider the following TITO process model of an industrial scale polymerization reactor developed by [33]: where the two outputs y 1 and y 2 are measurements representing the reactor conditions and the two inputs u 1 and u 2 are the set-points of two reactor feed flow loops. The proposed tuning method provides the PID controllers shown in Table 2. The upper bounds for this example are M S 1 = 1.60 for loop 1 and M S 2 = 1.20 for loop 2. These maximum sensitivity values are chosen to obtain a faster response for loop 1 and a smoother response for loop 2. The initial controller parameters are presented in Table 1.
For comparison, controller parameters are computed for the BLT [4] and Chen-Seborg [5] methods, and the results are also presented in Table 2. These methods are restricted to PI controllers. Additionally, a decoupler D = G −1 (0) must be considered for the Chen-Seborg case. Although the comparison is unbalanced when considering a decoupler for just one method, it shows the competitiveness of the proposed method even in this scenario.
The robustness and performance measures are provided in Table 3. The IAE values of each loop subjected to a set-point change are given in Table 4. As can be seen, the proposed method presents lower values of the IAE index compared to the other methods.
The Nyquist plots of the loop transfer functions are shown in Figure 4 for loop 1 and loop 2. Note that for the proposed method, the individual robustness index S ∞ is satisfied  for both loops. Moreover, the multivariable stability criterion is lower than that of the compared methods, which implies greater robustness when all loops are closed. The lowest cost criterion is achieved by the Chen-Seborg method, which is not surprising because of the use of the decoupler.
Unit set-point changes at time t = 0 h for y 1 and at time t = 10 h for y 2 are applied. The outputs y 1 and y 2 are shown  in Figure 5. For loop 1, the output response achieved by the proposed controller has less overshoot and a shorter settling time. The disturbance in loop 1 due to a set-point change in loop 2 is better rejected by the proposed method than by the BLT and Chen-Seborg methods. However, the disturbance rejection in loop 2 due to a set-point change in loop 1 is better for the Chen-Seborg controller; this is due to the fact that the Chen-Seborg case is the only case in which a decoupler is used. Figure 6 shows the controller outputs u 1 and u 2 due to unit step changes in the set-points. According to this figure, the Chen-Seborg method has the most aggressive control action (highest SII, see Table 4). Meanwhile, the proposed method has a smoother control response, requiring less from the actuator.
Note that although Figure 5 and the values of in Table 3 indicate superior performance of the Chen-Seborg method in terms of both output tracking and disturbance rejection due to the coupling effects, the authors would like to emphasize that a decoupler was used to achieve the calculated results. The fact that the achieved by the proposed method is 44 % greater than the of the Chen-Seborg method without a decoupler while, in similar conditions, the of the BLT method is 141 % greater than that of the Chen-Seborg method exemplifies how competitive the proposed technique is. VOLUME 9, 2021

B. EXAMPLE 2
In this example, the Shell heavy oil fractionator benchmark developed by the Shell Company [34] is considered. This multivariable process is highly constrained with strong interactions among its control loops and large dead times. In [35], the original system was modified by relaxing some of its constraints. Thus, it used the transfer function matrix of a 3 × 3 system given by Eq. (15). The three output variable are the top end pointy 1 , side draw end pointy 2 and bottom reflux temperaturey 3 ). Meanwhile, the three input variable represent the top drawu 1 , side drawu 2 and bottom reflux dutyu 3 . The physical details of this process can be found in [35].
The initial controller parameters are presented in Table 5. The chosen upper bounds for the robustness constraint are M S 1 = 2.30, M S 2 = 2.30, and M S 3 = 1.20. For loops 1 and 2 the maximum sensitivity values are chosen to obtain a faster response. Table 6 shows the controller parameters for the proposed method as well as the PI parameters of the method formulated by [35], which is referred to as the Lawal-Zhang method in this work.
The measures of robustness and performance are presented in Table 7. The performance cost function ( ) obtains higher values for the proposed method, 21.13 % higher compared to the Lawal-Zhang method result. For both methods, the multivariable stability criterion ( ) is fulfilled, i.e., ≤ 2n. The IAE and SII values of each loop are given in Table 8. Although the proposed controllers C 1 and C 3 have higher IAE values, note that the proposed method has lower SII values for all loops. These results indicate that the control action obtained with the proposed method is less aggressive than that with the Lawal-Zhang method.   The Nyquist plots of the loop gain transfer functions are shown in Figure 7 for loops 1 to 3. According to the Nyquist plots, the individual loop constraints are all satisfied for the proposed method. Note that the Lawal-Zhang method does not guarantee the stability of loop 1 individually. The stability of loop 1 is guaranteed only with loops 2 and 3 closed. This instability issue of loop 1 does not occur for the proposed method. From a practical point of view, this is an advantage of the proposed method. Figure 8 shows the time responses to unit step changes in the set-points. As can be seen, the proposed controller has a smoother response than the Lawal-Zhang method. Figure 9 shows the controller outputs due to unit step changes in the set-points. For all loops, the proposed method has a smoother control response (lower SII) compared to the Lawal-Zhang method.

C. EXAMPLE 3
Consider again the process given by (5). The Wood-Berry binary distillation column plant is a multivariable system that has been extensively studied. By applying the proposed method, two PID controllers were designed for the robustness constraints M S 1 = 1.70 and M S 2 = 1.70. These maximum sensitivity values are chosen to obtain a similar robustness condition for both loops. The initial controller parameters are presented in Table 9.
For comparison, the controller parameters were also calculated by the Boyd-Åström method described in [23]. The results of both methods are presented in Table 10.
The measures of robustness and performance are presented in Table 11. The individual loop constraints are all satisfied for the proposed method, which is also indicated in the Nyquist plots in Figure 10. The performance cost function of the proposed method is 42.25 % lower than that of the Boyd-Åström method. The log-modulus criterion ( ) results are similar for both methods. The performance indices of each loop subjected to a set-point change are given in Table 12. For both controllers, the proposed method results in lower values of the IAE and SII indices.     The closed-loop responses to unit step changes in the set-points for y 1 (at t = 0) and y 2 (at t = 80 min) are shown in Figure 11. The simulation results indicate that the set-point response in the y 1 output is faster for the proposed VOLUME 9, 2021  method than that for the PID controller tuned using the Boyd-Åström method. The disturbance response is significantly better for the proposed method when the set-point of loop 2 is changed. For the loop 2 output, the set-point and disturbance responses are similar for both methods. Figure 12 shows the controller outputs due to unit step changes in the set-points. The proposed method has a smoother control response (lower SII) compared to the Boyd-Åström method.

VII. CONCLUSION
In this paper, we have proposed a method for tuning decentralized PID controllers for multivariable systems. The tuned parameters of the decentralized PID controllers are obtained by solving a nonlinear optimization problem with two robustness constraints. The proposed method can be applied to higher-order multivariable processes that have complex dynamics and considerable multiple time delays. Furthermore, in contrast to the design and implementation difficulties of PI/PID controllers arising from decomposing MIMO systems into single loops, the proposed algorithm does not require such tedious work, and only one design parameter is necessary to specify the trade-off between robustness and performance. Three simulation examples are included to demonstrate the effectiveness of the proposed controller, and it can be observed that the aim of reducing the disturbance effects on the outputs due to set-point changes in coupled loops is accomplished.
Although the proposed technique yields strong robustness and stability results for the examples presented here, it still has room for improvement for industrial applications in the sense that it does not directly consider the effects on input and output constraints. In future works, we intend to extend the proposed technique by including these constraints in the optimization problem. This task is not trivial because a detailed study must be carried out to consider whether current constraints and additional process constraints overlap or make the optimization problem unfeasible. In addition, we intend to assess in which cases the designed controller parameter values are globally optimal.