Singularly Perturbation Method Applied To Multivariable PID Controller Design

1Department of Control and Mechatronics Engineering, Faculty of Electrical Engineering, Universiti Teknologi Malaysia, 81310 Skudai, Johor, Malaysia 2Department of Industrial Systems Engineering and Design, Jaume I University of Castello, 12080 Castello de la Plana, Spain 3Department of Industrial Electronics, Faculty of Electronics and Computer Engineering, Universiti Teknikal Malaysia Melaka, 76100 Durian Tunggal, Melaka, Malaysia


Introduction
Multivariable PID Control.Among the controller variety, PID becomes the controller that is most applied in a physical system [1].The reason is that it has a characteristic that offers simplicity, clear functionality, and ease of use [2].However, Ho et al. [3] reported that only one-fifth of PID control loops are in good condition.The others are not, where 30% of PID controllers are not able to perform well due to lack of tuning parameters, 30% due to the installation of a controller system operating manual, and 20% due to the use of default controller parameters.
In recent years, many researchers have paid attention to the MPID controller design for various systems such as in Industrial Scale-Polymerization Reactor [4], Coupled Pilot Plant Distillation Column [5], Narmada Main Canal [6], Quadruple-Tank Process [7], Boiler-Turbine Unit [8], and Wood-Berry Distillation Column [9].A research by Kumar et al. [4] had proposed a synthesis method of PI controllers based on approximation of relative gain array (RGA) concept to multivariable process.The method was further improved by relative normalize gain array concept (RNGA).Controller based on RNGA concept provides a better performance than RGA concept.Both concepts use the nonstandard PID controller which requires Maclaurin series expansion [10].In the work by Sarma and Chidambaram [5], PI/PID controllers based on Davison and Tanttu-Lieslehto method extended to nonsquare systems with right-half plane zero were applied.Results show that the Davison method gives better performance with less settling time than Tanttu-Lieslehto method.However, both methods are not applicable for square system.
Essentially, there are many integral controllers that are designed for nonlinear system [11,12].However, most of the existing techniques do not guarantee the desired transient performances in the presence of nonlinear parameter variations and unknown external disturbances [13].In a previous study, Martin and Katebi [14] had proposed Davison, Penttinen-Koivo, Maciejowski, and Combined method as a control tuning design for ship positioning.These controls strategies are based on PID controller which is used to control multivariable system.Due to the effectiveness and simplicity of those proposed controllers, Wahab et al. [15] had used those methods as tuning strategies for wastewater treatment plant (WWTP).The controllers were designed based on steady state of a linear system and static model inverse.The reliability of the proposed method was tested to a nonlinear WWTP.The response shows that good result was obtained from Davison until the Combined method.In the work by Balaguer et al. [16], a comparison between MPID controller with figure of merit controller was done for WWTP based on open-loop, closed-loop, and open-closed loop controller structure analysis.The MPID control tuning based on Davison and Penttinen-Koivo method was carried out by minimizing the residuals of both controllers obtained from the data.However, the dynamic nature of WWTP which involves ill condition characteristic causes difficulties in finding the optimum MPID tuning parameter.The system's behavior that involves slow and fast variables causes the control tuning strategies to not easily meet specification for multiple control variables at the same time.
A lot of approaches have been proposed to control multivariable system.Some of the approaches are able to deal with a high order multivariable system.However, a simple controller design has always become a desired controller where it is certainly can be accepted by the industry.By that, the required cost to run the system will be minimized as well.Realizing the simple controller design by other researchers [14][15][16], those methods were applied in this project and improved by adopting singularly perturbation method (SPM) to the controller designs by considering the dynamic matrix inverse.
Singularly Perturbed Multivariable Controller.Analysis and synthesis of singularly perturbed control have received much thoughtfulness over the past decades by many people from numerous arenas of studies [17][18][19][20][21]. Singularly perturbation method is able to decompose and simplify the higher order of the full order system into slow and fast subsystems [17,22,23], which are known as singularly perturbation system.Definitely, most of the control systems are dynamic, where the decomposition into stages is dictated by multitime scale.In this situation, the slow subsystem corresponds to the slowest phenomena and the fast subsystem corresponds to the fastest phenomena.It basically has two different parts of eigenvalue represented for slow and fast dynamic subsystems [24], where slow subsystem corresponds to small eigenvalue and fast subsystem corresponds to large eigenvalue [25].
This work is focused on the analysis of singularly perturbation system on two different case studies given.Singularly perturbed control of multivariable system is comprised of two steps.First, the multivariable system is decomposed into slow and fast subsystems.Then, the optimal composite singularly perturbed controller is designed [25][26][27][28].There are many approaches that have been developed concerning the control of singularly perturbation system.The approaches use different conditions on the properties of the used functions, different assumptions, different theorems, and different lemma [13,21,29] which are specifically based on the systems behaviour.
In a study by Rabah and Aldhaheri [24], singularly perturbation system has been modelled by using a real Schur form method.It shows that any two-time scale system can be altered into the singularly perturbed form via a transformation into an order real Schur form (ORSF).It is based on two steps, transformation of matrix A into an ORSF using an orthogonal matrix and then application of balancing algorithm to an ORSF.Li and Lin [17] had addressed the composite fuzzy multivariable controller to nonlinear singularly perturbation system.The composite controller was obtained from the combination of slow and fast subsystems.It was tested to a DC motor driven inverted pendulum system and it provides realistic and satisfactory simulation results.Multivariable control by Kim et al. [30] used successive Galerkin approximation (SGA) method.This method causes the complexity in computations to increase with respect to the order of the system.Therefore, singularly perturbation method was adopted to decompose the original system into slow and fast subsystems.Result shows that the use of the method greatly reduces the computation complexity and it is more effective than the original SGA method.
To the best of author knowledge, there are two methods to obtain the singularly perturbation system, which are by analytical [21,[29][30][31] and linear analysis [32][33][34][35].Singularly perturbation system obtained based on linear analysis is discussed and has been applied in this research.By exploiting the properties of singularly perturbation system to the dynamic matrix inverse of MPID control tuning methodology, an easy multivariable tuning method should be established.In Section 2, the time scale analysis is presented to determine the behavior of the system.Section 3 described the methods to obtain singularly perturbation system based on Naidu and Jian Niu method.The sequences of MPID controller based on Davison, Penttinen-Koivo, Maciejowski, and Combined methods are discussed in Section 4. Section 5 presented the optimization method which is based on particle swarm optimization (PSO).The case studies and the performance of the proposed methods for two case studies are investigated and discussed thoroughly in Sections 6 and 7. Finally, conclusions are given in Section 8.

Time Scale Analysis
To apply singularly perturbation method to the controller designs, the considered system must consist of a two-time scale characteristic.The two-time scale characteristic can be determine by rearranging the eigenvalue of the system in increasing order which will give  () = { 1 , . . .,   ,  1 , . . .,   } , (1a) where (), (  ), and (  ) are a total, slow, and fast eigenspectrum of the system, respectively. 1 is a smallest eigenvalue of the slow eigenspectrum,   is a largest eigenvalue of the slow eigenspectrum,  1 is a smallest eigenvalue of the fast eigenspectrum, and   is a largest eigenvalue of the fast eigenspectrum: The system is said to possess a two-time scale characteristic, if the largest absolute eigenvalue of the slow eigenspectrum is much smaller than the smallest absolute eigenvalue of the fast eigenspectrum.This is proven by where  is a measure of separation of time scales that represents an intrinsic property.

Singularly Perturbation Method (SPM) for MIMO System
Industrial processes possess "" number of inputs and outputs variables, where interaction phenomena exist.Interaction phenomena that occur among the inputs and outputs variables of multivariable process cause great difficulties in MPID controller design.Usually, it is solved by tuning the most important loop whereas other loops are detuned by keeping the interactions of that loop adequate.To compensate the interaction phenomena, one of the loops is forced to perform poorly.This detuning method is far from the optimal [4].In this work, to account for the interaction phenomena, instead of using the original process transfer function to the MPID controller design, that transfer function is rearranged by separating the slow and fast eigenvalues using SPM.

Naidu Method.
In this section we propose a procedure for a separation of slow and fast subsystem.The considered linear equations for two-time scale continuous system with the output vector possessing two widely separated groups of eigenvalues are where  and  are slow and fast variables in  and  dimension and  is a measured output.Matrices   ,   , and   are constant matrices of appropriate dimensions.Consider the problem as in (4a) to (4c).The system possesses a twotime scale property.Preliminary to separation of slow and fast subsystem, the system consists of  number of small eigenvalue (close to the origin) for slow subsystem and  number of fast eigenvalue (far from the origin) for the fast subsystem.The number of slow and fast eigenvalues needs to be identified based on eigenvalue location.Fast eigenvalue of the system is only essential during a short period of time.Then, it is insignificant and the characteristic of the system can be described by degenerating system known as a slow subsystem.
By letting  = 0, slow subsystem is obtained as By assuming  22 as a nonsingular matrix, (5b) becomes Using equation ( 6) in (5a), ẋ slow is represented as where Using ( 6) in (5c),  slow is represented as where To obtain fast subsystem, it can be assumed that the slow variables are constant at fast transients, where From (4b) and (6), Let

Mathematical Problems in Engineering
The fast subsystem is obtained as The composite system which consists of slow and fast subsystem is achieved using two-stage linear transformation which can be referred in an article written by Chang [36]: where The state space form of composite system is represented in (15a) to (15f).

Jian Niu Method.
The two-time scale system can also be solved using other method.This section presents singularly perturbation method based on Jian Niu.In order to apply Jian Niu method, transfer function matrix should be transform into a state space model: To illustrate the two-time scale decomposition, ( 16) is considered.Equation (16) can have this form where (  ,   ,   ) is a state space form of   ().Equations (17a) and (17b) can be represented as where  is a very small positive constant.Equations ( 18) can be considered as Equations (19a) to (19c) are the linear equations for two-time scale continuous system, similar just like (1a) to (1c) where This method is discussed in several literatures [33,35].The slow subsystem is denoted as And the fast subsystem is where The transfer functions for slow and fast subsystem are denoted by (30) and (31), respectively.The composite model is signified as a sum of slow and fast subsystem and a very little item () [37]  composite () =  slow () +  fast () −  slow +  () , (32) where

Multivariable PID Controller Design
Owing to the industrial process control involved with multivariable system, MPID controller design technique is necessary.It is a powerful control technique for solving coupling nonlinear system [38].The conventional MPID controller designs technique is based on static inverse model [15].This technique is difficult to obtain the desired control performance.Therefore, an enhancement is presented based on the dynamic inverse matrix and singularly perturbation method to the designs of MPID controller.Essentially, this enhancement has been considered in the previous work reported in [39] and it shows that the enhancement is able to control dynamic system where the output is able to track the set point change and produced less proses interaction.Nevertheless, it only considered that three controller designs instead of four and the selection of parameter tuning are done without optimization technique.In this paper, there are four enhanced MPID controller designs which are Davison, Penttinen-Koivo, Maciejowski, and Combined method where it is applied to the both original and singularly perturbed system with all of the parameter tuning being obtained based on particle swarm optimization.All of these four designs technique are applied to wastewater treatment plant and Newell and Lee evaporator.

Davison Method. Multivariable control design based on
Davison method simply applies the integral term, which causes decoupling rise at low frequencies The controller expression is represented by (34) [16], where   and () are integral feedback gain and controller error, respectively, Since this research is focused on dynamic control,   is defined as in (35), where  is the only controller tuning parameter, which undoubtedly needs to be tuned progressively until the finest solution is discovered.Due to the involvement of the inverse system, the control design is only applicable for square matrix.If the system involves time delay, the time delay needs to be eliminated.

Penttinen-Koivo
Method.This method is somewhat advanced and then the method proposed by Davison.In Penttinen-Koivo method, a proportional term is introduced.It comprises both integral and proportional term.Indirectly, this causes decoupling to take place at low and high frequencies.Davison and Penttinen-Koivo method are only similar in terms of an integral term which is linearly related to the inverse of plant dynamics The controller expression is represented in (36) [16], where   ,   , and () are proportional gain, integral feedback gain, and controller error correspondingly Dynamic terms of   and   are expressed in (37), where  and  are the tuning parameters for both proportional and integral feedback gain.

Maciejowski Method.
Maciejowski method applies all proportional, integral, and derivative gains in its controller design.For maciejowski method, the tuning was done around the bandwidth frequency,   .Consequently, the interaction is reduced and good decoupling characteristic is provided around the frequency [15].However, due to the needs of plant frequency analysis experiment, this method is quite difficult to be used throughout the industry [15] The controller expression is represented by (38), where   ,   ,   , and () are proportional, integral feedback, derivative gains, and controller error Dynamic terms of   ,   , and   are expressed in (39), where , , and  are Maciejowski tuning parameters.Due to a complex gain obtained from the calculation of (  ) −1 , a real approximation of (  ) −1 is necessary which can be done by solving the following optimization problem: where  is a constant that is used to minimize .

Combined Method.
In order to overcome the weakness of the Maciejowski method which requires rigorous frequency analysis, a new method was proposed by Wahab et al. [15].It is the result of the previous controllers where methods by Davison, Penttinen-Koivo, and Maciejowski are combined together: Equation ( 41) represents the proposed control design, where , , and () are the tuning parameters and controller error: is defined in (42). is also a tuning parameter.This method keeps some properties in Maciejowski method but excludes the needs of frequency analysis [15].

Optimized Singularly Perturbed MPID Parameter Tuning
To ensure a fair comparison, the optimum parameter tuning for each of controller designs is measured by using particle swarm optimization (PSO).PSO optimizes a problem by having a population (swarm) of candidate solutions (birds) which is known as particles that are updated from iteration to iteration [40].These particles are moved into the search space seeking for a food according to its own flying experience and its companion flying experience.It can be expended to multidimensional search.Each particle (solution) is characterized by its position and velocity, and every one of them searches for better positions within the search space by changing its velocity [41].Each particle preserves the track of its current position within the search space.This value is identified as the particle's local best known position (pbest) and leads to the best known position (gbest), which is defined as enhanced positions that are found by the other particles.By that, the finest solution is attained Equation ( 43) represents the update equations of new velocity and new position, where V( + 1), ( + 1), , V(),  1 ,  The overall performance of PSO can be increased by proper selection of inertia weight, .Lower value of  provides a good ability for local search and higher value of  provides a good ability for global search [41]: To achieve a respectable performance,  is determined according to (44), where  max is the maximum value of inertia weight,  min is the minimum value of inertia weight, iter is the current number of iteration, and iter max is the maximum number of iteration.Most of the previous researchers have used  max = 0.9 and  min = 0.4, where significant enhancement of PSO is achieved [42,43]: Fitness function which is also known as cost function is represented by (45), where () is a system error.The procedure of PID parameter optimization by using PSO is summarized as follows: (1) Initialization: initialize a population of particles with arbitrary positions and velocities on X dimensions in the problem space.Then, randomly initialize pbest and gbest.(2) Fitness: calculate the desired optimization fitness function in X dimensions for every particle.(3) pbest: compare calculated fitness function value for every particle in the population.If current value is smaller than pbest, and then update pbest as current particle position.(4) gbest: determine the best success particle position among all of the individual best positions and designate as a gbest.
(5) New velocity and position: update the velocity and position of the particle based on (43).
(6) Repeat step (2) all over again until a criterion is encountered.

Case Studies
In the next subsections, an introduction to the case studies is presented.First, an overview of the wastewater treatment plant (WWTP) is provided and the Newell and Lee evaporator model is explained.These two case studies are considered to demonstrate the performance of the proposed methods.

Case Study I: Wastewater Treatment Plant (WWTP).
Wastewater treatment plant is designed either for carbon removal or for carbon and nitrogen removal.In this project, the carbon removal scenario is considered.The control plant outputs are substrate and dissolved oxygen.Scanty provision of substrate affects the growth of microorganisms that are responsible for treating the wastewater and too many provisions lead to a drop in the microorganisms growth rate.The standard amount of substrate is around 51 mg/L [44].Meanwhile, insufficient dissolved oxygen will cause the degradation of the pollutants and the plant to become less efficient.Too much dissolved oxygen can cause excessive consumption of energy where it will increase the cost for the treatment.Other than that, it also can cause too much sludge production.The amount of dissolved oxygen concentration needs to be controlled so that it is in the range of 1.5 mg/L-4.0mg/L [45].The aim of this case study is to control the concentration of substrate and dissolved oxygen at the desired value by manipulating the manipulate variables of dilution rate and air flow rate, respectively.The state space form of the nonlinear wastewater treatment plant is linearized from [39] as follows: where the state is composed by , the biomass, , the substrate, DO, the dissolved oxygen, and   , the recycled biomass.
These , , , and  matrices are used in the design of singularly perturbed MPID control and tested into the nonlinear wastewater treatment plant.

Case Study II: Newell and Lee
Evaporator.This subsection presents the Newell and Lee evaporator system which is considered as a second case study.The objective is to evaluate the effectiveness and the performance of the proposed singularly perturbed MPID controllers for different system.Here, unstable system is considered.Similar to the first case study, the four different methods of MPID are implemented, which is Davison, Penttinen-Koivo, Maciejowski, and Combined methods.The plant to be controlled is given by the following state space model [46]: where the state is composed by  2 , the separator level,  2 , the product composition, and  2 , an operating pressure.The input variables are  2 , the product flow rate,  100 , the steam pressure, and  200 , the cooling water flow rate.The outputs to be controlled are

Results and Discussion
This section presents the results and discussion for both case studies.It will include the results of singularly perturbed MPID control based on four proposed methods, which are Davison, Penttinen-Koivo, Maciejowski, and Combined method for both full order and singularly perturbed system.The first section shows the MPID control based on particle swarm optimization.The obtained optimum tuning parameters are presented.The second section shows the simulation results for control of the closed loop system during substrate and dissolved oxygen set point change, while the last section provides the result which shows the stability of the closed loop system.

Results for the
By using algorithm discussed in Section 3.1, the original system which refers to (47a) to (47d) is represented in state space form of singularly perturbed system as indicated in (52a) to (52d).This state space form is used to design the controller tuning, while the simulation and performance of the system are based on the original system.From the state space form in (52a) to (52d), it is clearly shown that the eigenvalues of the system are grouped into two distinct and separate sets, which causes the time consumed to obtain the MPID tuning parameters reduce.All eigenvalues of the singularly perturbed system are remained at the left-half plane, which is −1.9955, −0.2798, −0.0214, and −0.0660.It is indicated that the system is boundary input boundary output (BIBO) stable Based on the algorithms explained in Section 3.2, the Jian Niu method is successfully able to represent the original system into composite of singularly perturbed system as in (53a) to (53d), where the only existing eigenvalue is located at −0.0508.Since it is in the left-half plane, the system is BIBO stable.To validate the models from both methods, singularly perturbation method obtained by Naidu and Jian Niu, the magnitude and phase plot between the original system and singularly perturbed system are plotted as shown in Figure 1.
Based on Figure 1, singularly perturbed system based on Naidu algorithm provides better dynamic response compared to the singularly perturbed system based on Jian Niu algorithm, where the dynamic response is much identical to the original system.The close approximation between the original system and singularly perturbed system by Naidu in the frequency responses analysis exhibits the authority of model and essentially leads to adequate control performance of the controller design.Hence, a singularly perturbed system by Naidu is used in dynamic matrix inverse of MPID controller designs.
In the considered WWTP, there are two controlled variables and two manipulated variables.Interaction phenomena may occur between these two controls and manipulate variables.Each manipulated variable can affect both the control variables.The process interaction among the variables may cause the closed loop system to become destabilized and the controller tuning is more challenging.In order to minimize the process interaction, the selection of suitable control and manipulated variables pairing is importance.In this case, there are two possible controller pairings.
The relative gain array (RGA) analysis has been used in quantifying the level of interactions in a multivariable system.It is also used to determine the best input output pairing and that pairing should be avoided.To measure the ability of RGA in providing a realistic pairing recommendation, the RGA for To out(2) To out (2) From in (1) From in(2) From in(1) From in(2) linearized models was calculated.RGA of a nonsquare matrix is defined in (54).The results are displayed in a matrix form, where columns are for each input variable and rows for each output variable.This matrix form can be used in determining which relative gains are associated to which input output variables where From the RGA obtained, it can be concluded that dissolved oxygen cannot be paired with dilution rate due to the negative relative gain.It corresponds to the worst case, and this is highly undesirable.Biomass, substrate, and recycle biomass on the other hand cannot be paired with air flow rate due to the zero relative gain, which means that air flow rate does not have any effect on biomass, substrate, and recycle biomass.
In RGA analysis, the closer the value of RGA element to one is, the configuration is more likely to work, where less interaction exists.Hence, it is concluded that a good pair of dissolved oxygen and air flow rate and substrate and dilution rate are recommended.Since MPID controller designs are involved with several tuning parameter, PSO was adopted.Due to the PSO characteristic which cannot give a unique solution at every attempt [47], 10 trials of simulation for original and singularly perturbed system of each MPID controller design were conducted.A result with minimum error was selected.Table 1 shows the obtained optimum PID tuning parameter using ITSE fitness function for both systems: original and singularly perturbed system.The results are corresponding to Davison, Penttinen-Koivo, Maciejowski, and Combined method, respectively.Based on the results presented, it clearly shows that singularly perturbed system is able to provide easiness in tuning strategy in terms of computation time where it required less computation time compared to the original system.Table 1 shows that singularly perturbed MPID based on each method is able to reduce more than half of computation time required by original system.
Figures 2 to 5 show the comparison between output response and interaction based on Penttinen-Koivo method for original and singularly perturbed system.Based on the Figures 2 and 4, Penttinen-Koivo based on singularly perturbed system is able to produce better output responses  with less overshoot and fast settling time.Meanwhile, more oscillation is produced by output responses based on original system.The interaction among the output variables for Penttinen-Koivo based on singularly perturbed system is also reduced as shown in Figures 3 and 5.
Due to the better closed loop performance, reduced process interaction, and less time consuming obtained by the Penttinen-Koivo method based on singularly perturbed system compared to the original system, singularly perturbed system was implemented thoroughly in the case of the controller performance evaluations among others three methods accordingly.To measure the performance quality of four singularly perturbed MPID controller designs, pseudorandom binary sequence (PRBS) was injected as the input signal to the system.PRBS was injected to determine and to test the tracking ability of the proposed singularly perturbed MPID controller designs for each step change.From the responses obtained, all four designs are able to track each step change where the Combined method provides the best response.Figures 6 and 7 show the output responses for both substrate and dissolved oxygen concentration.
To provide a clear view of the set point tracking ability and process interaction, the system was also injected with     step input.The simulation was carried out during substrate and dissolved oxygen set point change.For each change, the step input was injected at  = 10 h and  = 50 h, respectively.Figures 8 and 9 show the simulation results for substrate and dissolved oxygen for each proposed singularly perturbed MPID controller design.The responses are set with respect to the step change in the substrate input from 41.2348 mg/L to 51.2348 mg/L.
Figure 8 shows that all singularly perturbed MPID controller designs are able to keep the concentration of the substrate close to the desired value.It shows that the control based on Combined method is able to provide the finest control effect among the others' method in terms of  settling time and maximum amplitude, where it is able to achieve settling point during 11 h compared to the others which settle during 44 h, 21 h, and 14 h, respectively.Due to the control characteristic which only applies integral gain, control action based on Davison method provides a response with the highest percentage of overshoot (%OS).By using Penttinen-Koivo method, the output response shows better improvement.The presence of both integral and proportional gain is able to minimize the percentage of overshoot (%OS) and offer a better settling time (  ).However, proportional gain needs to be tuned wisely.High value of proportional gain can cause the system to become unstable, while small value of proportional gain may reduce the sensitivity of the controller.Compared to the Davison and Penttinen-Koivo method, Maciejowski method gives the best performance with small percentage of overshoot (%OS) and faster settling time (  ).It is proven that control performance at the selected frequency was able to improve the closed loop response.An important feature in Maciejowski method is the selection of frequency.Frequency must be selected properly to avoid instability.Among all methods, the Combined method exhibited the best tracking to the substrate changes.This method exhibits a faster response than other control designs, but it requires a long time to obtain the tuning parameters.
Since the considered case study involves multivariable system, process interaction may occur.Interactions between the system variables occur because each manipulated variable in the multivariable system certainly will affect the controlled variables.Here, dilution rate and air flow rate will affect the response of both substrate and dissolved oxygen concentration.Evidently, when changing one of the inputs for dilution rate or air flow rate, both outputs will be affected, and this means that there is significant coupling in the system.Figure 9 shows the interaction responses for each proposed singularly perturbed MPID controller design that occurs during substrate change.The response obtained has proven that substrate and dissolved oxygen are coupled since the step changes in the substrate disturb the dissolved oxygen correspondingly.If there is no process interaction, dissolved oxygen should not be affected when the substrate is changed.Fortunately, process interaction was reduced for each controller design where the Combine method provides less interaction, which indicates by the lowest maximum amplitude and less oscillation.Large variation exhibits in the response of dilution rate and air flow rate based on Davison method.Penttinen-Koivo and Combined method exhibit moderate variations, while Combined method exhibits smallest variations but with high peak value at one of the points.
Figures 12 and 13 show the simulation results for dissolved oxygen and substrate to the step change in the dissolved oxygen input.From Figure 12, it is clearly shown that the trajectory of the output is improved by each method that has been applied.It is proportionate with the responses during substrate set point change.All the systems step response settles at a final value of 4.1146 mg/L, which is the final value of the unit step input.Figure 12 shows that the responses produced by the Davison and Penttinen-Koivo methods do not asymptotically approach the final value, where overshoot appears in the final value.The maximum values of the outputs are 4.8649 mg/L and 4.2830 mg/L, respectively, for each of Davison and Penttinen-Koivo methods.The output response yield by Davison method consists of 37.52% of overshoots.However, Penttinen-Koivo method has improved the output response by reducing the overshoot to 8.42%.Meanwhile, the responses acquired from Maciejowski and Combined methods asymptotically approach the final value.These methods provide slightly similar effect in terms of maximum amplitude and settling time.There is no overshoot of the final value and there are no oscillations in the response.The outputs reach the final value at around  = 0.2 h and  = 0.4 h for both Maciejowski and Combined methods, respectively.Figure 13 shows that interactions also occur during dissolved oxygen set point change.Similar to the responses during substrate set point change, process interaction exists and it is improved by each method proposed by Davison up to the Combined method, respectively.Figures 14 and 15 show the closed loop responses of manipulate variable, which are dilution rate and air flow rate during the dissolved oxygen set point change, respectively.The characteristic of the closed loop responses for each of singularly perturbed MPID controller design is summarized in Table 2.
The stability of a system can be determined directly from its transfer function or from CLCP. Figure 16 shows the closed loop poles and zeros plot for each singularly perturbed MPID control design.It is mark a pole location by a cross (x) and a zero location by a circle (o).Based on the plot figure, all poles are located on the left-half plane that guarantees a stable system.However, to ensure the reliability of the stability analysis, Routh-Hurwitz analysis was performed.The results show that all methods are able to produce a stable system.As a result Since  is less than 1, the system behaved as a two-time scale characteristic.There is one slow variable which is indicated by eigenvalue of 0 and two fast variables which are indicated by the eigenvalue at −0.0558 and −0.1000.Based on the algorithm discussed in Section 3.1, the original system of Newell and Lee evaporator can be represented in singularly perturbed system.The eigenvalues for singularly perturbed system are 0, −0.0558, and −0.1000, which is similar to the original system From in (1) From in (2) From in (1) From in (2) To out (1) To out (2) To out (1) To out(2) Figure 17: Bode analysis for original and singularly perturbed system.
Equations (59a) to (59d) represent the singularly perturbed system of Newell and Lee evaporator in state space form.
To verify the singularly perturbed system with the original system, the Bode diagram is plotted as shown in Figure 17.
As seen from Figure 17, singularly perturbed system shows a fairly good tracking with the original system over the low, middle, and high frequencies ranges.The response of the singularly perturbed system is almost identical to the original system.The close approximation between both systems demonstrates the validity of the obtained singularly perturbed system and principally leads to satisfactory control performance.
In this case, there also exists two possible control and manipulate variables paring.By using (54), the RGA for Newell and Lee evaporator was obtained as Based on the analysis of RGA, it can be concluded that separator level cannot be paired with cooling water flow rate, and operating pressure cannot be paired with product flow rate, respectively.This is due to the zero relative gain.A change of cooling water flow rate will not give any significance to the separator level, and change of product flow rate will not give any significance to the operating pressure.From the RGA analysis, it is highly recommended to pair product flow rate with separator level and cooling water flow rate is paired with operating pressure.Since the value is nonzero and positive, the pairing is possible.Similar as in case study I, 10 trials of PSO simulation for original and singularly perturbed system of each MPID controller design were conducted.However, due to the unstable open loop response, Maciejowski method cannot be implemented to this multivariable system as it requires information from stable open loop response.Table 3 shows the obtained optimum PID parameter based on PSO for evaporator system.By applied Davison method, both systems are able to provide similar tuning parameter with similar error for 10 number of run.However, it can be seen the advantage of singularly perturbed system which required less computation time compared to the original system.For Penttinen-Koivo and Combined method, the computation time is reduced more than triple times with the adaptation of singularly perturbed system in MPID control.Figures 18 to 21 show the comparison between output responses based on Penttinen-Koivo method for each Newell and Lee evaporator, original and singularly perturbed system.Figure 18 shows the separator level responses between original and singularly perturbed system during separator level change.The response belonging to the original system has a poor performance as compared to the singularly perturbed system with high oscillation during step up and step down response.The output response achieved by the singularly perturbed system is with less overshoot and fast settling time.Figure 19 shows the interaction response during separator level change.It can be seen that the interaction was reduced by the adaptation of singularly perturbed system in MPID control.Figure 20 shows the operating pressure response between original and singularly perturbed system during operating pressure change.It can be observed that the singularly perturbed system has better performance as compared to the original system with fast settling time, while Figure 21 shows the interaction response during operating pressure change.Due to the good responses exhibited from the  control based on singularly perturbed system, this system was implemented thoroughly into the Davison, Penttinen-Koivo, and Combined method.The simulation was performed by setting the set point input values as follows: value of separator level: 1 m, value of operating pressure: 50.5 kPa.
The simulation was executed in two-variable change, which is during separator level and operating pressure step point change.The step point change of the separator level and operating pressure is listed in Table 4.
Figures 22 and 23 show the closed loop performance for the separator level and operating pressure to the sequential step point changes in separator level set point.The step point changes were sequentially introduced into the system at  = 120 s and  = 250 s, respectively.In the simulation study, the comparisons of the closed loop performances were done between  which is based on multiobjective optimization approach using surrogate modelling.Fauzi's method was compared as they also have been involved with the similar Newell and Lee evaporator system.However, the details concerning to the controller designed are not presented and can be referred in [46].Figure 22 shows that the output performance based on Davison method provides a response with 11.67% and 8.02% overshoot during the step up and step down input, respectively.These values are higher compared to other controller methods.However, it is still able to track the set point input given.The output performances based on Penttinen-Koivo and Combined methods are quiet similar.Combined method provides faster rise and settling time during the step up input, whereas the Penttinen-Koivo provides faster rise and settling time during the step down input.However, the Combined method provides the best performance with low percentage of overshoot, which specify by the lowest value of maximum amplitude and steady state error.Among the four methods, the method proposed by Fauzi is the poorest.At the early stage, the response is relatively good.Once the step down input is injected, the response shows unstable characteristic where it fails to settle at the given set point input.The response is considered as unstable since the gain error increases as time increases.Figure 23    large interactions with high maximum amplitude and more oscillation.The sluggishness in the performance is due to the controller algorithm which only involves integral gain.It can be observed that the Combined method is able to reduce the interaction effects well compared to the other methods.The interaction produced by the Combined method is the lowest.Penttinen-Koivo method produces interaction slightly higher than the Combined method, while Fauzi method offers quiet high interaction.
Figures 24 and 25 show the closed loop responses of manipulate variable, which are product flow rate and cooling water flow rate during the separator level set point change respectively.For Davison method, large variations of product flow rate and cooling water flow rate are obtained, while Penttinen-Koivo and Combined method consist of small variations but high peak value.Among the four methods, Mohd Fauzi method exhibits the largest variations.slightly high steady state errors than the Combined method.The difference is only about 1.40%.Meanwhile, the difference is approximately 0.1% during step down input.Even though the response by Davison method required long computation time for the rise and settling, the response is accomplished to settle at the set point value.But the response is relatively slow and consists of high percentage overshoot and steady state error.Figure 26 also shows that the output performance based on Davison method provides a response with 10% overshoot during the step up and step down input.These values are higher compared to other controller methods.However, it is still able to track the set point.The output performances based on Penttinen-Koivo and Combined methods are almost similar.Penttinen-Koivo provides faster rise and settling time during the step up and step down input.However, the Combined method provides the best performance with the lowest steady state error.Meanwhile, the resulting response by Fauzi obviously shows unstable characteristic.At the beginning, the response already shows that the system is in a state of uncontrollable.After the step down input was injected, the response gradually decreased.At time  = 800 s, the response of operating pressure is at −366.7 kPa.An increase in simulation time will lead the response to be infinity.Figure 27 shows the response of interactions during operating pressure set point change.Among the four methods, Penttinen-Koivo and Combined methods offer the least interaction.It can be seen that the interaction is reduced with the Penttinen-Koivo and Combined methods compared to the Davison and Fauzi methods which consist of high maximum amplitude.
Figures 28 and 29 show the closed loop responses of corresponding manipulated variable, which are product flow rate and cooling water flow rate during operating pressure set point change, respectively.The variation of both manipulate variables is similar during separator level set point change, where control based on Mohd Fauzi method exhibits a response with the largest variations.
The characteristic of closed loop response for evaporator system of all the comparative methods for singularly perturbed MPID control during the separator level and operating pressure set point change is tabulated in Table 5.The good performance of the Combined method is readily apparent.The best performance is given by Combined method, followed by Penttinen-Koivo, Davison, and Fauzi method.
Figure 30 shows the closed loop pole-zero plots for the proposed singularly perturbed MPID controller designs applied to the Newell and Lee evaporator.It can be seen that  all eigenvalues are located at the left-half plane of the s-plane.This indicates that the closed loop system is generally stable.A real pole in the left-half plane defines an exponentially decaying component in the homogenous response.The rate of the decay is determined by the eigenvalue location.Eigenvalues far from the origin in the left-half plane correspond to the components that decay rapidly, while eigenvalues near the origin correspond to slowly decaying components.Referring to Table 5, the rise and settling time of the response based on Penttinen-Koivo method are the most faster.It is proportional to the poles location indicated in Figure 30.

Conclusion
Designing MPID control tuning based on original and singularly perturbed system for multiinput multioutput (MIMO) processes is presented.Simulation results lead to the inference that, with the appropriate parameter tuning, a satisfactory singularly perturbed MPID control performance can be accomplished to control a nonlinear model of wastewater treatment plant.Ill-defined system like wastewater treatment plant which usually faces difficulties in control system, due to the natural behavior of two-time scale characteristic, can be efficaciously controlled by the implementation of singularly perturbed system into the MPID controller designs.Among the four methods, the Combined method yields somewhat better results with respect to decoupling capabilities, closed loop performances, and process inteaction.
For the second case study, Davison, Penttinen-Koivo, and Combined method were successfully applied to the nonlinear model of Newell and Lee evaporator.The well-tuned parameters of the controller designs were obtained using PSO approach.Simulation results show that the implementation of singularly perturbed system to the dynamic matrix inverse of Davison, Penttinen-Koivo, and Combined method has

Figure 1 :
Figure 1: Bode analysis for different methods.

Figure 2 :
Figure 2: Substrate responses between original and singularly perturbed system during set point change.

Figure 3 :
Figure 3: Process interaction between original and singularly perturbed system during substrate change.

Figure 4 :
Figure 4: Dissolved oxygen responses between original and singularly perturbed system during set point change.

Figure 5 :Figure 6 :
Figure 5: Process interaction between original and singularly perturbed system during dissolved oxygen change.

Figure 7 :Figure 8 :
Figure 7: Concentration of dissolved oxygen based on PRBS input.

Figure 9 :
Figure 9: Process interactions during substrate set point change.
Figures 10 and 11  show the closed loop responses of manipulate variable, which are dilution rate and air flow rate during the substrate set point change, respectively.

Figure 10 :Figure 11 :
Figure 10: Dilution rate responses during substrate set point change.

Figure 12 :
Figure 12: Closed loop responses of dissolved oxygen during dissolved oxygen set point change.

Figure 13 :Figure 14 :
Figure 13: Process interactions during dissolved oxygen set point change.

Figure 15 :Figure 16 :
Figure 15: Air flow rate responses during dissolved oxygen set point change.

Figure 18 :
Figure18: Separator level responses between original and singularly perturbed system.

Figure 19 :
Figure 19: Process interaction between original and singularly perturbed system during separator level change.

Figure 20 :
Figure 20: Operating pressure responses between original and singularly perturbed system.

Figure 21 :
Figure 21: Process interaction between original and singularly perturbed system during operating pressure change.

Figure 22 :Figure 23 :
Figure 22: Responses of separator level during separator level set point change.

Figure 24 :Figure 25 :
Figure 24: Product flow rate responses during separator level set point change.

Figure 26 :
Figure 26: Responses of operating pressure during operating pressure set point change.

Figure 27 :
Figure 27: Process interactions during operating set point change.

Figure 28 :
Figure 28: Product flow rate responses during operating pressure set point change.
The input variables are , the dilution rate, and , an air flow rate.Matrices , , , and  are given by

Table 1 :
Optimum PID parameter for WWTP based on PSO.

Table 3 :
Optimum PID parameter for evaporator system based on PSO.

Table 4 :
Step point change.
Davison, Penttinen-Koivo, Combined, and multivariable controller designed proposed by Fauzi, compares the interactions that occur during the separator level set point change.It is clearly shown that Davison method produces