PLC based fractional-order PID temperature control in pipeline: design procedure and experimental evaluation

Fractional-order control system design can be used for systems with non-local dynamics involving long-term memory effects. However, implementation of a fractional-order controller in industrial systems is complicated, because of the excessive demand for computational power. The following paper presents the step-by-step design procedure, parameter tuning, and experimental evaluation of the fractional order proportional-integral-derivative (FOPID) controller. The control algorithm is based on the Continued Fraction Expansion approximation of the fractional-order operators. It is implemented on a standard industrial Programmable Logic Controller. The FOPID control system is verified and evaluated in terms of efficiency and robustness using a new laboratory benchmark of a temperature control in the pipeline. The proposed solution shows increased efficiency in terms of robustness compared to the standard PID closed-loop control.


Introduction
The Proportional-Integral-Derivative (PID) algorithm is one of the most used algorithms in the industry to control processes. It has a simple implementation and there is a vast number of easy-to-use and intuitive tuning methods for such a controller depending on the controlled process with required control quality [20,38,48].
Many physical phenomena exhibit complex microscopic behaviours. Most of the processes associated with complex systems have non-local dynamics involving long-term memory effects. To model such processes with higher accuracy, fractional calculus, i.e. the integration and differentiation of the arbitrary order, can be used [28,32,39]. Fractional derivative operators are non-local. This means that a functional derivative of function f(x, t) depends on all the values of this function in the interval ½t 0 ; t (entire history of f(x, t)). Some examples of using fractional calculus for modelling physical processes are e.g. acoustics, diffusive transport, fluid flow, geodynamics, optics, chemical engineering and electromagnetic waves [16].
The standard PID control algorithm and dynamic models used in control engineering are based on integer order calculus (a special case of fractionalorder calculus) [8]. It can be shown that integer derivative operators don't hold the non-locality property [11]. To effectively control the dynamic systems described by fractional-order equations, fractionalorder control algorithms should be designed. The major problem with control engineering applications of fractional calculus was the complicated method of solving the differential equations. It was overcome with finding the proper geometric and physical interpretation of the fractional operators [56] and their efficient numerical approximations, like Power Series Expansion (PSE) [54] or Continued Fraction Expansion (CFE) [7], among others. One family of such controllers is the set of Fractional Order PID-type (FOPID) controllers, such as PID k [31], CRONE (fr. Contrôle Robuste d'Ordre Non Entrier, Non-Integer Order Robust Control) [47], and PI k D l [55]. Such controllers have increased the number of tuned parameters from 3 to 5 (compared with standard PID) and allow better matching of the controller with the process. There are also different tuning methods proposed for the FOPID controller which fulfill the control system design requirements e.g. optimization with integral criteria [9,32], constrained min-max optimization [5], swarm optimization [21], autotuning methods [33], and robust tuning methods [59]. In [9,32,65] authors show that the FOPID controller outperforms the standard PID controller in applications where processes exhibit non-local dynamics.
In the literature, many practical applications of the FOPID controller and it's variants are proposed. In [10] authors present FOPI controller for wind turbine generators. Applications of FOPID control in electrohydraulic systems are considered in [18,60,66]. In [3] the twin-rotor helicopter control system based on FOPID controller is described. The optimal tilt FOPID control in rail vehicles is presented in [15]. In [27,35] the implementation of FOPI controller in industrial electric drives (DC-motors and PMSM drives) is described. Application of FOPID control in the pumpturbine governing system (PTGS) is proposed in [26]. FOPID control of the Automatic Voltage Regulator (AVR) is given in [50,64]. In [61] FOPID control system for pumped storage hydro unit is described.
FOPID fuzzy control of the combined cycle power plant is proposed in [14]. In [25] the FOPID controller is implemented in the pneumatic servo-system. In [58] the FOPID controller is applied in Load Frequency Control (LFC) in power systems. In [29] a decoupled visual model and a FOPID controller are designed aiming at the problem of the view distortion in the fillet weld welding process.
Most FOPID algorithms are usually tested and simulated in a MATLAB-Simulink environment [10,14,15,26,50,58,61,62,64]. Controllers that implement the fractional-order control algorithms usually are equipped with fast and efficient DSP processors (e.g. dSpace 1 ), because extensive computations are required [18,27,29,32,35,66]. It makes such a solution hard to apply in the industry. The tuning of controller parameters is also a problem that needs to be solved to ensure the proper functionality of the automation system. Nowadays, most production lines in small and medium companies are controlled by Programmable Logic Controllers (PLC). Therefore, PLC implementation is favorable in industrial applications. In [23] authors present the PLC implementation of a CRONE 3rd generation controller for a hydraulic system based on SIEMENS S7-300 PLC. The PLC implementation of the FOPID algorithm with CFE approximation for temperature control using a BR 2005 controller is described in [51]. The FOPID algorithm with PSE approximation for the water level in a hydraulic system is given in [36] . Implementations of fractional operators based on CFE approximation using SIEMENS S7-1200 PLC are presented in [42], and Hardware-in-the-Loop (HIL) experiments are described in [43]. Finally, preliminary results of laboratory benchmark for robust PLC-based FOPID temperature control in the pipeline is presented in [53].
Many solutions presented in the literature on the effective PLC FOPID implementation, are tested in the Hardware-in-the-Loop setup [36,[42][43][44][45]. Only few describe the laboratory experiments with real processes [25,36,51]. This article synthesizes the mentioned works and describes the step-by-step design procedure, parameter tuning, and implementation of the PI k D l algorithm with CFE approximation on a laboratory test stand. The control system is equipped with a Siemens S7-1200 PLC for control and data acquisition. The controlled value is the temperature in the pipeline with induced air-flow. Moreover, on the test stand, it is possible to impose additional disturbances, to check the robustness of the proposed solution. In the paper the comprehensive tests, covering both the set-point tracking and the disturbance rejection, are considered, which is rarely mentioned in the articles about industrial PLC implementations of FOPID controllers. Therefore, robustness of the PI k D l algorithm can be thoroughly tested. The proposed solution shows increased efficiency in terms of robustness compared to the standard PID control.
The article is organized as follows: In Sect. 2, basic fractional-order calculus definitions are given. Then, in Sect. 3 the PI k D l controller algorithm and its parameters' tuning are described. In Sect. 4 numerical approximation of fractional operators is explained. Next, in Sect. 5, the laboratory test stand is described. In Sect. 6 the the PLC implementation of the PI k D l controller is presented. In Sect. 7 experimental results for the PI k D l controller, along with a comparison with the standard PID controller are shown and discussed. Finally, concluding remarks are given in Sect. 8.

Differ-integral operator
The fractional calculus generalizes the well-known integer order actions by using real numbers for differentiation and integration actions. The differintegral operator is the combination of non-integral differentiation and integration operators, defined as where a-order of the operation (a 2 R), a and t are the lower and upper limits of integration. The differ-integral operator (1) acts as the integrator for a\0 and as the derivative for a [ 0.
The Riemann-Liouville definition (RL definition) of a differ-integration of function f(t) is given as [11] RL a D a t f ðtÞ ¼ d a dt a f ðtÞ where m À 1 a m, m 2 N-Riemann-Liouville fractional derivative-integral, CðÁÞ-Euler Gamma function [19], defined as The most common differ-integral operator definitions are the Riemann-Liouville (RL definition) [11], the Caputo (C definition) [6], and the Grünwald-Letnikov definition (GL definition) [28,32].
Because of the need of numerical implementation in PLC controllers, usually the Grünwald-Letnikov definition of a differ-integration of the function f(t) is given as where a j ¼ is a generalization of the Newton symbol into real numbers. For numerical simulations, the GL definition (4) can be rewritten as GL a D a t f ðtÞ ¼ where w ðaÞ j -binominal coefficients, calculated recursively as follows It should be emphasized, that RL (2) and GL (4) formulations are equivalent [30]. The Laplace transform of (2) can be calculated as [54] L RL 0 D a t f ðtÞ where FðsÞ-Laplace transform of the function f(t), s-complex variable, m [ 0; m À 1 a\m

Oustaloup approximation of differ-integral operator
Direct numerical implementation of the differ-integral operator (4) is difficult to apply, and it is necessary to approximate it. The commonly used approximation method in practical implementations is the Oustaloup filter in the form of a product of a series of stable firstorder linear systems [4,46] s The Ostaloup filter coefficients are calculated as follows where x d , x g -lower and higher bound frequencies, n-the order of the approximation.
3 PI k D l controller algorithm and tuning The proposed PI k D l controller transfer function is in the following form [55] where K P FOPID -the proportional gain, T I FOPID -the integral time constant, T D FOPID -the derivative time constant, s k , s l ,-differ-integral operators (k\0, l [ 0). The block diagram of the controller structure is given in Fig. 1, where the input of the controller is the difference between the measured process value and set-point value given as where y m ðtÞ-the measured plant output, y SP ðtÞ-setpoint. PI k D l controller can be tuned according to the optimization method described in [49,52,65], assuming the controller transfer function in the form (14) and the process transfer function in the form The chosen tuning method is based on the Oustaloup approximation (9). It requires information about the desired gain margin A m and phase margin h m of the closed-loop control system to solve the following set of equations where x g , x p -the gain and phase crossover frequencies.
To calculate x p , x g the following condition has to be fulfilled During the controller parameters calculations, to fulfill the conditions (21)-(23), the numerical optimization methods are usually used. Then, if the parameters x p , x g , k, l are given, the controller gains K P FOPID , T I FOPID , T D FOPID can be calculated from (17)-(20) as follows 4 Discrete approximation of differ-integral To implement the differ-integral operator in the digital controller (e.g. PLC), different discrete approximation methods are used. The first common approximation method is the Power Series Expansion (PSE) [42] that implements directly GL differ-integral. PSE leads to a solution in the form of polynomials, that have a digital Finite Impulse Response (FIR) filter structure. Such an approximation of the differ-integral operator depends only on the zeroes of the filter, so only system outputs are considered in the solution. Therefore, reasonable PSE approximation requires output signal values from many previous steps, which increases the requirements for the controller's hardware and memory. Another approximation method is the Continued Fraction Expansion (CFE) [40,42]. CFE approximation leads to the solution in the form of rational functions. Therefore, it has an Infinite Impulse Response (IIR) filter structure. It depends on both zeros and poles of the system, so it requires previous outputs and inputs of the system. With CFE, there is a lower requirement for output and input signal values from previous steps required than with PSE. It allows for fast convergence and low order of the filter.
For the discrete approximation of the fractional operator, the generating function is used in the form where z-complex variable. Then, the CFE approximation of the fractional operator can be written as follows The discrete transfer function, approximating the fractional operator is given as where a 2 ½0; 1-argument that depends on the approximation method (a ¼ 0 for Euler approximation and a ¼ 1 for Tustin approximation), T p -sampling period, Y(z)-discrete transform of the output sequence yðnT p Þ, F(z)-discrete transform of the input sequence f ðnT p Þ, p, q-orders of the approximation (usual assumption is p ¼ q ¼ n, so then n is the approximation order). The time response of the element (29) is expressed as where y(k), u(k)-the output and input signals in provided sequence, v i and w j -coefficients of CFE approximation (for exact values see [40]), k-discrete Since CFE provides better memory usage, it can be used for approximation of fractional operators in real systems with limited capacity.

Description of laboratory test stand
The process part of the laboratory test stand is presented in Fig. 2, whereas the schematic diagram is given in Fig. 3.
The laboratory test stand comprises a pipeline with a heater (G) in which airflow is induced by a fan (S). During experiments, the control system changes the amount of heat released by the heater (G) and maintains constant air-flow. The power control signal of the heater (Y G ) and speed of the fan (Y W ) are standard current signals 4-20 mA generated by a PLC controller used in a closed-loop control system.
To measure the temperature a T/I measuring transducer with a Pt100 resistance thermometer and additional linearization was used. The measurement range of the transducer was 25 C -75 C (output standard current signal 4-20 mA). To control the system, a PLC S7-1212c Siemens controller was used with an analog input/output module and a 24V power supply. A SIMATIC HMI KPT 600 panel and desktop PC with TIA-Portal software were used for visualization and controller implementation.
It was possible to artificially induce the following disturbances: • a step-change in the air inlet orifice -change of the cross-section from 389 to 1661 mm 2 (DI1); • a step-change in the heater power (G) by adding or disconnecting additional resistance which changes the heater resistance from 100 X to 75 X (DI2). 6 PLC Implementation of PI k D l controller The design and PLC implementation procedure of the PI k D l (FOPID) controller is presented in Fig. 4. It consists of 7 steps described below.

Process identification
The mathematical model of the process was approximated by the transfer function of the First Order Plus Dead Time (FOPDT) system.
where PVðsÞ ¼ yðsÞ-the process value, CVðsÞ ¼ uðsÞ-the control value, T z -the process time constant, T 0 -the process dead time, K ob -the process gain.
The plant parameters were determined from the open-loop system step response using the secant method [12]. The FOPDT model parameters were calculated as follows where DCV-step change of CV, DPV-amplitude of the step response, T 1 , T 2 -time constants for the step reponse to reach 50% and 63,2% of its final value. Then, the Second Order model of the process was calculated. It allows for getting the transfer function in the form (16), based on the FOPDT transfer function (31). The Second Order model parameters can be approximated as follows [22] where The goodness of fit measure between the estimated output signal y m ðkÞ and the reference output signal y r ðkÞ was calculated as Normalized Root Mean Square Error [17] NRMSE ¼ 1 À ky est ðkÞ À y r ðkÞk ky r ðkÞ À y r k ; ð35Þ where y est ðkÞ-estimated output signal, y r ðkÞ-reference signal, y r -mean value of the reference signal, k-sample. Value of NRSME varies between À1 (bad fit) to 1 (perfect fit).
Both identified models, with corresponding NRMSE values, are gathered in Table 1. In Fig. 5 step responses are presented and in Fig. 6 estimation errors for identified models are presented. For both models quality index is around 0.98, which indicates a high modeling accuracy.
The Interior-Point optimization method [37] with the Integral Squared Time Error (ITSE) criterion was used to fulfill the conditions (21)- (23). The calculated parameters of the PI k D l controller, according to (24)- (26), are gathered in Table 2.

Controller implementation
The implementation process of the PI k D l control algorithm (14) on PLC can be divided into two phases.

Phase 1: simulink implementation
In the first phase, the PI k D l control algorithm is implemented in the MATLAB/Simulink environment. The set of MATLAB functions and functional Simulink blocks should be created. In Fig. 7 the overall Simulink diagram of the implemented PI k D l controller is given.
In the block CFE (see Fig. 8) the Continued Fraction Expansion approximations of the fractional operators s k and s l are implemented. In blocks I_fopid and D_fopid (see Figs. 9, 10) the  Fig. 5 Step responses comparison for the identified models   Fig. 7) realizes the controller transfer function (14).

Phase 2: transfer to PLC
In the second phase, the code for the PLC is implemented in the Structured Control Language (SCL). For this purpose, Simulink PLC Coder 2 can be efficiently used, which allows code creation for PLC directly from the MATLAB code. The program in the TIA Portal environment (Siemens) can be divided into blocks. In this way, a modular structure of the application is obtained. There are 4 types of blocks: • Organizational Blocks (OB) allow the code to run when a specific condition occurs. The condition may be, for example, a measured period, a specific date, error, or hardware configuration change. • Function Blocks (FB) the instructions they contain are also executed cyclically. However, the instructions they contain usually need over one program loop cycle. Therefore, it is often necessary to know the function values from the previous cycles. For this purpose, data blocks (DB) are used. • Functions (FC) used to perform simple and repetitive instructions such as assignment instructions or simple arithmetic calculations, • Data Blocks (DB) these blocks make it possible to store data determined during controller operation, and are also available in every other block.
Implemented TIA program structure, with Function Blocks realizing FOPID algorithm and accompanying Data Blocks, is presented in Fig. 11. The controller code has been divided into 5 function blocks. The fbCFE determines the coefficients according to the CFE approximation described. It takes the sampling period, fractional integration factor k, fractional derivative factor l, and the factor determining the discretization method as the parameters. As a result of this function, 4 tables are created that are the weighting factors of the subsequent inputs  Finally, the controller response is calculated in block fbFOPID. For the program to work, it was also necessary to implement an organizational block that starts at the time of starting the controller to determine the initial conditions and the block calling cyclic interrupts with a frequency equal to the assumed sampling frequency.

Quality indices
During tests of the controller, the quality of control has to be evaluated. It was analyzed in the time domain using the following indices [34] • Steady-state error e st -error in steady-state after introducing the set-point change or the disturbance. • Maximum error where T 2 ½T sp ; T R -it is the time-span between the T sp -first moment of reaching the set point value and settling time where e 1 and e 2 are the first two consecutive biggest errors with opposite signs, assuming as the zero level (baseline) the steady-state value of the output signal y(t) after the transient response.
• Settling time T R -it is the time between the moment of change of the set point y sp ðtÞ, or introduction of disturbances d(t) and the moment when error e sp ðkÞ reaches a fixed value inside a boundary e.g. AE0:1%.

Selection of the sampling period
Digital devices such as a PLC or a micro-controller have to work with a finite sampling period (T p ). Thus, the process output is available only at the sample points and the process variable is discrete. To guarantee the expected performance in terms of the continuous plant output it is crucial to choose an adequate value of T p . When discretizing continuous controllers by using mathematical formulas describing the relation between Laplace and Z-transforms, a smaller sampling period means higher accuracy. Numerous research was done in this topic (considering constant and variable sampling rates) [24,41,63]. For the FOPDT system (31), the following practical rule is proposed [2] T p % T z N r ð38Þ where N r -a number of sampling periods during T z (rise time). In PLC implementations, the value of N r is usually chosen from the range 4-10 [2]. The equation (38) allows calculating the maximum sampling period. For the process considered in this work T z ¼ 43:6691 s (see Table 1), thus the sampling period should be lower than 4 s (N r ¼ 10). There is no theoretical lower limit for the sampling period for the considered controller. However, PLC with a too-short T p cannot provide all required calculations in real-time, and better performance, though it can lead to unnecessary costs of high-end instrumentation and capacity resources. On the other hand, a low sampling rate will reduce the control quality because the controller's response to interference is too slow.
We have chosen the sampling period T p ¼ 200 ms for PID and FOPID for simulations and PLC implementation. It should be mentioned, that sampling in the Siemens S7-1200 PLC is not fixed (as in MATLAB/Simulink), and there exist AE0:005 s margin. This T p value was chosen concerning the requirements for computing performance (with the I/O manipulation, HMI communication or data monitoring and acquisition) and control quality. There was observed, that increasing further the sampling rate decreased quality because some of the control cycles didn't produce proper control value.

Simulations
During simulations, the algorithm implemented in the PLC was used to control the process simulated in MATLAB/Simulink environment (Software-In-the-Loop, SIL). For this purpose, it was necessary to transfer data between the Siemens PLC software (TIA Portal) and the MATLAB/Simulink environment. An additional OPC UA server was used to enable data exchange between the simulated controller in the PLCSim (part of the TIA-Portal software) and Simulink. NetToPLCsim software 3 extension allowed TCP / IP connections of PLC with the simulated process. Calculated quality indices for SIL tests were gathered in Table 3. This confirmed the effective transfer of code from the MATLAB/Simulink environment to the Siemens S7-1200 PLC controller.

Laboratory tests
During laboratory tests, the designed and tuned PI k D l controller was used in the laboratory test stand. The step response of the control system with the PI k D l controller, and the influence of disturbances were examined.
For the comparison, the integer-order PID controller was chosen, in the form where K P PID -proportional gain, T I PID -integral time constant, T D PID -derivative time constant and The integer-order PID controller was tuned with a Ziegler-Nichols method [67], which is well known among control engineers. These tuning rules give significant disturbance rejection. However, the drawback is the possibility to observe some overshoot in the step response. Calculated PID controller settings are gathered in Table 4.
All experiments started with the operating point: The calculated control quality indices for the control system with the PI k D l controller are summarized in Table 5.

Discussion
The results shown in Tables 3 and 5 prove successful implementation of the PI k D l algorithm in the PLC controller. The proposed system evaluated experimentally on a laboratory test stand, allowed tracking of the set-point with high accuracy. It was also robust to induced disturbances (air inlet change, and heat power change).
The first major problem was a PI k D l tuning method that requires off-line optimization, and more parameters to choose than a standard PID controller. The tuning rules are based on the Oustaloup continuous approximation.
Because of implementation requirements in the PLC controller, the discrete CFE approximation of fractional operators was used. Therefore, there still exists the problem of interpreting optimized controller parameters for a PLC-based discrete PI k D l control system. When compared to the SIL simulation, the laboratory implementation was characterized by a four times shorter settling time and no overshoot. It was a surprisingly good result. However, the experimental comparison with the PID controller allows observing other important system properties that are hard to detect in simulations.   In all evaluated cases, the maximum errors didn't exceed 0.4% (approx. 0:2 C) for PI k D l , and 1.25% (approx. 0:6 C) for PID controller. For the set-point change, an important advantage of the PI k D l controller was zero overshoot and shorter settling time compared with the PID controller. As it was mentioned before, the occurrence of overshoot for the PID controller depends on the Ziegler-Nichols tuning method properties.
Both PID and PI k D l controllers, after proposed parameter tuning, are robust to disturbances. In this case, the advantage of the PI k D l controller is around two times smaller maximum error comparing to the standard PID controller. However, it has a worse settling time and higher overshoot. Such differences in quality indices for a PI k D l controller can be described  ). This phenomenon is well known for controllers with derivative action, which amplifies abrupt changes in the error signal (input to the controller). For the PID controller, it also occurred (see Figs. 18 and 22), but the magnitude of noise was smaller. There are many ways to improve signal quality in closed-loop control systems. Methods like noise filtering, feedforward from load disturbances, or process dynamics compensation are used in industrial practice [13,57]. However, there is no systematic research on such methods for fractional-order control systems.

Conclusions
The article describes the design, parameter tuning, and experimental evaluation of a fractional-order PI k D l temperature control in a pipeline with induced air-flow on a laboratory test stand. The PI k D l controller was implemented on a standard PLC Siemens S7-1200 controller. Fractional operators were approximated using the Continued Fraction Expansion (CFE). The proposed solution was tested via simulations and experiments on the laboratory test stand, and compared with a standard integer-order PID controller.
The results obtained during the presented research, with the thermal-flow process, showed that to apply the proposed control algorithm it is crucial to develop methods for tuning the controller parameters to minimize the influence of input delays and disturbances. The proposed future work includes research on tuning methods of discrete FOPID controllers with CFE approximation of differ-integral and improvements of signal quality in closed-loop fractional-order control systems.

Compliance with ethical standards
Conflict of interest The authors declare that they have no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.