# Identification of Continuous-Time Systems with Multiple Unknown Time Delays from Sampled Data

Authors: Yamna Ghoul; Naoufel Zitouni; Kaouther Ibn Taarit; Moufida Ksouri
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

The purpose of this work is to identify a multiple-input single-output (MISO) continuous time (CT) system with various unknown time delays using sampled data. A novel technique for simultaneously estimating system parameters and various time delays may be developed using a gradient algorithm. Indeed, we propose a new formulation of the identification problem which permits defining the generalized vector containing the time delays and the dynamic parameters. The gradient optimization algorithm with exact derivatives of the objective function with respect to the adjustable parameters has been suggested as an alternative for solving this optimization problem. The identification method can be easily implemented by a filtering procedure to generate the input/output time derivatives. In this sense, the proposed method, which is based on a state variable filter, intrinsically allows data filtering and thus simplifies the use of the method in a practical case. An analysis is then presented which proves the convergence of the algorithm, and robustness issues are discussed as well. Finally, the suggested technique is evaluated on simulated instances to ensure its validity and performance.


---

## Body

## 1. Introduction

Time delays are unavoidable in the operation of many industrial process models [1, 2]. Indeed, ignoring its presence might result in complicated behavior (oscillations, disregarded data, instability of the closed loop, and so on) [3, 4]. As a result, the identification of time delay models is indeed a problem of considerable importance that deserves special attention in virtually all disciplines of science and engineering. Two major issues must be addressed: the first step is to identify the dynamics parameters, and the second step is to identify the temporal delays. Indeed, many techniques take into account both the identification of CT systems [5, 6] and the identification of DT systems [7, 8].In fact, the identification of physical systems was usually conducted in discrete time. This was mainly due to the “go completely discrete-time” trend that was spurred by parallel developments in digital computers.Indeed, historically discrete-time models have received more attention than continuous-time models. Therefore, many researchers in the system identification field have concentrated their efforts on developing DT identification methods [7, 9–11]. In fact, for some reasons, discrete-time approaches have been recognized. First, the stochastic theory of DT approaches is more easily compared with CT approaches which need to calculate time derivatives of input/output data. Second, estimating the physical parameters of the continuous-time physical model is not important if the model is to design controllers.However, continuous-time approaches are preferable under certain circumstances. First, when we have to estimate physical signal process inputs such as pressure and temperature, we should use continuous-time approaches because the physical signal is always a continuous-time signal rather than a discrete-time signal from a zero-order holder. Second, when we need to estimate physical parameters for the purposes of analysis and system design, continuous-time approaches are recommended.Indeed, many continuous-time system identification methods dealing with CT systems with time delays have been developed [5, 12–15].In reality, the challenge of identifying SISO time-delay systems is addressed by the majority of these approaches. Nonetheless, one of the most challenging issues is identifying MISO time-delay systems that is a field of study where few works have been published. Nevertheless, for multi-input single-output (MISO) systems with multiple unknown time delays, the problem is much more difficult. Such systems are often more challenging to model. In particular, systems with several inputs can be difficult. A basic reason for the difficulties is that the couplings between several inputs and an output lead to more complex models. A second reason for studying the MIS0 systems with multiple unknown time delays is the presence of multiple unknown delays that is generally considered as a great challenge where few works have been developed. This decision is also supported by the fact that these models are often employed in most process sectors, including vehicle traffic [16], thermal agitation [17], network management [18], and variable pitch milling [19]. In summary, it is worth noting that the identification issue of MISO CT systems with multiple unknown time delays is, therefore, both of theoretical and practical importance.However, several studies have previously been conducted on the issue. In [20], a technique based on a modulating function method was proposed. Paplinski employs an ant colony optimization technique in [21]. This method, in reality, has a considerable execution time. The evolutionary algorithm [22], which estimates the optimal parameters, is another approach for estimating systems with many unknown time delays. This scheme is said to be an offline technique that takes too long to compute at times. Another technique is described in [23], which permits the detection of MISO CT systems with unknown time delays. There are two versions: one is based on least squares, while the other is not. The second version, on the other hand, is based on the instrumental variable estimator, which is utilized to provide unbiased estimates. GSEPNLS and GSEPNIV are the acronyms for global separable nonlinear least-squares and global separable nonlinear instrumental variable, respectively. In [24], the problem of MISO CT model identification with multiple unknown time delays from sampled data has been treated. The proposed method estimated the plant and the time delays in a separable way. More precisely, the plant has been estimated by the standard recursive least-square algorithm while the time delays have been explicitly estimated by the Gauss–Newton algorithm. In [25], we have estimated the parameters applying a Levenberg–Marquardt algorithm while the time delays have been estimated by the Gauss–Newton algorithm. In [26, 27], the separable estimation of continuous-time hybrid “Box–Jenkins” systems with multiple unknown time delays has been considered. It can be seen that the developed methods [24–27] require a separable identification strategy of parameters and time delays. This strategy is considered to be very slow and sometimes requires a large computational effort. However, successful implementations of identification methods in terms of computation time require simultaneous identification of the unknown time delays and the linear parameters of the system. An interesting issue, in this regard, includes simultaneously identifying the model parameters and time delays.The majority of the developed approaches require a separable identification strategy of parameters and time delays. This technique is focused on using a specific algorithm to identify temporal delays. The linear parameters are then calculated using a different approach. In truth, this method is quite sluggish and frequently necessitates a significant amount of computing work. However, in order to reduce the computation time, successful identification methods require simultaneous identification of the unknown time delays and the linear parameters of the system. An interesting issue, in this regard, includes simultaneously identifying the model parameters and time delays.In order to overcome the underlying problems, this study addresses the topic of simultaneous identification of the MISO CT system with various unknown time delays to address the underlying issue. Indeed, we provide a novel formulation of the issue that permits the temporal delays and dynamic parameters to be defined in the same estimated vector. Then, based on this formulation, we offer a novel gradient algorithm-based approach. The developed algorithm’s convergence criteria are also developed.The remainder of the present paper is laid out as follows. The problem in this study is introduced in Section2. The developed gradient algorithm is shown in Section 3. In Section 4, the suggested algorithm’s convergence is tested. We provide two instances in Section 5 to demonstrate the efficacy of the proposed approach. Finally, in Section 6, final observations are offered to bring the study to a close.


## 2. Problem Statement

In this research work, the following process is considered:(1)xut=G1sFse−d1su1t+⋯+GrsFse−drsurt.e−djs is the time delay operator defined as follows:(2)e−djsujt=ujt−dj.So,(3)xut=G1sFsu1t−d1+⋯+GrsFsudt−dr.The polynomialsGjs and Fs are defined, respectively, by(4)Gjs=gj0+gj1s+⋯+gjqjsqj,Fs=1+f1s+⋯+fpsp.The considered system can be described as(5)∑i=0pfisp−ixut=∑j=1r∑i=1qjgjisqj−iujt−dj.To sum up, the considered system is depicted in Figure1.Figure 1
The considered process.


## 3. Suggested Technique

Let us introduce the low-pass filterK(s) [28] where α and p are the cutoff frequency and the order.(6)Ks=1αs+1p.Using the bilinear transformation described in [29], we obtain(7)x¯0k+∑i=1pfiξix¯k=∑j=1r∑i=1qjgjiξp−qj+iu¯jk−d˜j+νk.Such that(8)νk=∑i=0pfiξiν¯k,ξiν¯k=ν¯kKiq−1,and(9)ξiu¯jk=u¯jkKiq−1,ξix¯k=y¯kKiq−1,Kiq−1=T/2i1+q−1i1−q−1p−iα1−q−1+T/21+q−1p.


### 3.1. Definition

d˜j is computed as(10)d˜j=djT=kj+ΔjT.The approximated DT identification model described in (7) is rewritten as follows:(11)x¯0k=φTk,d^θ^+vk,where θ^ and d^ are described as(12)θ^T=f^T,g^1T,...,g^rT,with(13)f^T=f^1,...,f^p,g^jT=g^j1,...,g^jqj,d^T=d^1,...,d^r.Consider the following prediction error denoted by(14)ek=x¯0k−θ^φTk,d^.φk,d^ is defined as(15)φTk,d^=−φx¯Tkφu¯1Tk−d˜1φu¯rTk−d˜r⋮φu¯rTk−d˜r,with(16)φx¯Tk=ξ1x¯k,...,ξnx¯k,φu¯jTk−d˜j=ξp−qj+1u¯jk−d˜j,...,ξpu¯jk−d˜j.We consider (17) which describe the generalized vector:(17)Θ^T=θ^T,d^TT.The following criteria is adopted:(18)VΘ^=1N−hs∑k=hs+1N12e2k.(N-hs) is the data length such that hs≥maxd˜1,...,d˜r.The generalized observation vector is the following:(19)Φk,Θ^=−∂ek∂Θ^=φk,d^,−∂ek∂d^.The gradient of the error is computed as in [30]:(20)∂εk∂d^j=−∑i=1qjg^ji∂∂dju¯jp−qj+ik−d^j|t=kT=−∑i=1qjg^jiddtu¯jp−qj+ik−d^j|t=kT=−∑i=1qjg^jiAdq−1T/2p−qj+i−1×1−q−1p−qj+i−11−q−1pj+i−1u¯jp−qj+ik−d^j=−∑i=1qjg^jiu¯jp−qj+i−1k−d^j.Finally, we obtain(21)Φk,Θ^=−∂ek∂Θ^=−φx¯Tkφu¯1Tk−d˜1⋮φu¯rTk−d˜r−∑i=1q1g^1iu¯jp−qj+i−1k−d˜j⋮−∑i=1qrg^riu¯jp−qr+i−1k−d˜r.The following gradient algorithm is used:(22)Θ^i+1=Θ^i−γiV′Θ^i.So,(23)V′Θ^=∂VΘ^∂Θ^=∂1/N−hs∑k=hs+1N1/2e2k∂Θ^=1N−hs∑k=hs+1N∂ek∂Θ^ek.(24) can be derived from (23):(24)V′Θ^i=−1N−hs∑k=hs+1NekΦk,Θ^.


## 3.1. Definition

d˜j is computed as(10)d˜j=djT=kj+ΔjT.The approximated DT identification model described in (7) is rewritten as follows:(11)x¯0k=φTk,d^θ^+vk,where θ^ and d^ are described as(12)θ^T=f^T,g^1T,...,g^rT,with(13)f^T=f^1,...,f^p,g^jT=g^j1,...,g^jqj,d^T=d^1,...,d^r.Consider the following prediction error denoted by(14)ek=x¯0k−θ^φTk,d^.φk,d^ is defined as(15)φTk,d^=−φx¯Tkφu¯1Tk−d˜1φu¯rTk−d˜r⋮φu¯rTk−d˜r,with(16)φx¯Tk=ξ1x¯k,...,ξnx¯k,φu¯jTk−d˜j=ξp−qj+1u¯jk−d˜j,...,ξpu¯jk−d˜j.We consider (17) which describe the generalized vector:(17)Θ^T=θ^T,d^TT.The following criteria is adopted:(18)VΘ^=1N−hs∑k=hs+1N12e2k.(N-hs) is the data length such that hs≥maxd˜1,...,d˜r.The generalized observation vector is the following:(19)Φk,Θ^=−∂ek∂Θ^=φk,d^,−∂ek∂d^.The gradient of the error is computed as in [30]:(20)∂εk∂d^j=−∑i=1qjg^ji∂∂dju¯jp−qj+ik−d^j|t=kT=−∑i=1qjg^jiddtu¯jp−qj+ik−d^j|t=kT=−∑i=1qjg^jiAdq−1T/2p−qj+i−1×1−q−1p−qj+i−11−q−1pj+i−1u¯jp−qj+ik−d^j=−∑i=1qjg^jiu¯jp−qj+i−1k−d^j.Finally, we obtain(21)Φk,Θ^=−∂ek∂Θ^=−φx¯Tkφu¯1Tk−d˜1⋮φu¯rTk−d˜r−∑i=1q1g^1iu¯jp−qj+i−1k−d˜j⋮−∑i=1qrg^riu¯jp−qr+i−1k−d˜r.The following gradient algorithm is used:(22)Θ^i+1=Θ^i−γiV′Θ^i.So,(23)V′Θ^=∂VΘ^∂Θ^=∂1/N−hs∑k=hs+1N1/2e2k∂Θ^=1N−hs∑k=hs+1N∂ek∂Θ^ek.(24) can be derived from (23):(24)V′Θ^i=−1N−hs∑k=hs+1NekΦk,Θ^.


## 4. Summary of Proposed Algorithm

The above approach can be summarized by the following step-by-step procedure:Step 1: initialization(1)
Leti = 0(2)
Set the initial value of the tuning parametersγ0, and the initial value of the cost function V0Step 2: gradient computationCompute the gradientV′Θ^i using (24)Step 3: gradient optimization algorithmComputeΘ^i+1=Θ^i−γiV′Θ^iStep 4: stopping criteriaIfΘ^i+1 has converged or the maximum number of iterations is reached, then stop, else increase i by 1 such that i = i + 1 and go back to Step 2End


## 5. Convergence Analysis

Theorem 1.
The following condition should be verified:(25)0<γi<2N−hs×1∑h=hs+1NΦk,Θ^ΦTk,Θ^.Proof.
ReplacingV′Θ^i in (22) with its expression (29), we obtain(26)Θ^i+1=Θ^i+γiN−hs∑k=hs+1NekΦk,Θ^.
Such that from (14), the prediction error is given as(27)ek=x¯0k−ΦTk,Θ^Θ^.
By replacinge(k) in (26) with its expression given in (27), we lead to(28)Θ^i+1=Θ^i+γiN−hs∑k=hs+1NΦk,Θ^x¯0k−ΦTk,Θ^Θ^.
Let denotesΘ^opt the optimal value of Θ^. We can then rewrite (32) as(29)Θ^opi+1=Θ^opi+γiN−hs∑k=hs+1NΦk,Θ^x¯0k−ΦTk,Θ^.
Or at the convergence, we have(30)Θ^i+1=Θ^i.
Consequently, we obtain(31)γiN−hs∑k=hs+1Nx¯0k−ΦTk,Θ^Θ^optΦk,Θ^=0.
It derives then(32)x¯0k−ΦTk,Θ^Θ^optΦk,Θ^=0.
Then,(33)Φk,Θ^x¯0k=Φk,Θ^ΦTk,Θ^Θ^opt.
So, we get(34)Θ^i+1=Θ^i+γiN−hs∑k=hs+1NΦk,Θ^ΦTk,Θ^Θ^opt−Θ^i.
Now, let define(35)χi=Θ^opt−Θ^i.
Then,(36)Θ^i+1=Θ^i+γiN−hs∑k=hs+1NΦk,Θ^ΦTk,Θ^χi.
We subtract the termΘ^opt from (36):(37)Θ^i+1−Θ^opt=Θ^i−Θ^opt+γiN−hs∑k=hs+1NΦk,Θ^ΦTk,Θ^χi.
Taking into account the expression ofχi in (35), it follows:(38)−χi+1=−χi+γiN−hs∑k=hs+1NΦk,Θ^ΦTk,Θ^χi.
Multiplying (38) by (−1), we obtain(39)χi+1=χi−γiN−hs∑k=hs+1NΦk,Θ^ΦTk,Θ^χi.
We get(40)χi+1=χi1−γiN−hs∑k=hs+1NΦk,Θ^ΦTk,Θ^.
The condition described in (41) should be verified to guarantee the convergence:(41)1−γiN−hs∑k=hs+1NΦk,Θ^ΦTk,Θ^2<1.
We obtain, finally,(42)0<γi<2N−hs×1∑k=hs+1NΦk,Θ^ΦTk,Θ^.
End of proof.Theorem 2.
The following theorem is verified:(43)EΘ^−ΘΘ^−ΘT=1N−hs2∑k=hs+1NΦk,Θ^ΦTk,Θ^δ2.Proof.
By the use of the first-order Taylor series expansion around the real parameter ofΘi, we obtain(44)V′Θ^i=V′Θi+V″Θi×Θ^i−Θi.
The proposed algorithm estimatesΘi satisfying V′Θ^i=0, so it can be derived from (44) that(45)Θ^i−ΘiΘ^i−ΘiT=V″Θi−1V′ΘiV′ΘiTV″Θi−1T.
It follows:(46)EΘ^i−ΘiΘ^i−ΘiT=EV″Θi−1V′ΘiV′ΘiTV″Θi−1T.
Or:(47)V″Θi=eke″k−e′kΦk,Θ.
Then,(48)V″Θi=eke″k+ΦTk,ΘΦk,Θ.
Using the small residual algorithms, the termeke″k can be neglected:(49)V″Θi≈ΦTk,ΘΦk,Θ.
The mean value ofV′ΘV′ΘT permits to get(50)EV′ΘV′ΘT=1N−hs2E∑k=hs+1NΦk,Θ^ΦTk,Θ^∑k=hs+1NekeTk.
Or,(51)E∑k=hs+1NekeTk=δ2.
Finally, (52) can be simplified as follows:(52)EΘ^−ΘΘ^−ΘT=1N−hs2∑k=hs+1NΦk,Θ^ΦTk,Θ^δ2.


## 6. Numerical Example

Two examples are presented to confirm the performance of the developed scheme.The vector of parameters to be estimated is the following:(53)Θ=f1,f2,g11,g12,g21,g22,d1,d2.


### 6.1. Example 1

The first simulated example is described as(54)x¨t+f1x˙t+f2xt=g11u˙1t−d1+g12u1t−d1+g21u˙2t−d2+g22u2t−d2.In the Laplace domain, the above equation becomes(55)xs=e−d1ss+1u1s+2e−d2ss+2u2s=e−d1ss+2u1s+e−d2s2s+2u2ss2+3s+2.The real parameters are the following:(56)f1=3,f2=2,g11=1,g12=2,g21=2,g22=2.d1=8.83s,d2=2.32s.


#### 6.1.1. Simulation Results

We obtain the following results illustrated in Table1.Table 1
Simulation results with noise-free output.
TrueEstimatesf132.9818 ± 0.0182f221.9844 ± 0.0156g1110.9912 ± 0.0088g1221.9876 ± 0.0124g2121.9941 ± 0.0059g2221.9830 ± 0.017d18.838.8301 ± 0.0001d22.322.3255 ± 0.0055With aSNR=20dB, the parameters estimated are shown in Table 2.Table 2
Results simulation withSNR=20dB.
TrueEstimatesf132.9813 ± 0.0187f221.9840 ± 0.016g1110.9931 ± 0.0069g1221.9887 ± 0.0113g2121.9955 ± 0.0045g2221.9810 ± 0.019d18.838.8276 ± 0.0024d22.322.3240 ± 0.0032


#### 6.1.2. Model Validation

Validations of the obtained model by using the gradient algorithm are plotted, respectively, in the following Figures2 and 3.Figure 2
Model validation of the proposed scheme with noise-free output.Figure 3
Model validation in the case ofSNR=20dB.


#### 6.1.3. Bode Diagrams

The bode diagrams are plotted, respectively, in Figures4 and 5.Figure 4
Bode diagrams of the true system model with the estimated model: noise-free case.Figure 5
Bode diagrams of the true system model with the estimated modelSNR=20dB.The accuracy of the developed prediction error method is gradually improved in the frequency domain.


#### 6.1.4. Computational Time

A comparison is done in Table3 with the computational time in the second simulation case. For comparison, the identification algorithm GSEPNIV for the global separable nonlinear instrumental variable given in reference [23] is also adopted to estimate the above stable-type model parameters defined by (53) by taking the same initial parameter setting.Table 3
Computational timeSNR=20dB.
GSEPNIVProposed methodComputational time (sec)5.7307581.756012The computational time is measured using the tic and toc Matlab function.Table3 shows that the GSEPNIV method needs a long computational time. However, it can be seen that the computational time is smaller by using the proposed prediction error method.


#### 6.1.5. Noise Effects

Table4 addresses another comparison between the proposed method and the GSEPNIV method developed in [23] in terms of robustness with regard to noise in the case of a noisy output.Table 4
Parameters estimates: simulation 2.
TrueGSEPNIVProposed methodf133.0646 ± 0.06462.9813 ± 0.0187f222.0773 ± 0.07731.9840 ± 0.016g1110.9804 ± 0.01960.9931 ± 0.0069g1222.0759 ± 0.07591.9887 ± 0.0113g2122.0269 ± 0.02691.9955 ± 0.0045g2222.0758 ± 0.07581.9810 ± 0.019d18.848.8671 ± 0.02718.8276 ± 0.0024d22.322.3315 ± 0.01152.3240 ± 0.0032Judging from the simulation results tabulated in Table4 concerning computational error values for each parameter, it can be seen that the proposed method outperforms the GSEPNIV method.


### 6.2. Example 2

We consider the following process described in Figure6, where u1 is the steam flow, u2 is the water flow u2, and x is the water temperature.Figure 6
The steam-water heat exchanger.The equation of the above process is(57)x¨t+f1x˙t+f2xt=g11u˙1t−d1+g12u1t−d1+g21u˙2t−d2+g22u2t−d2.Then,(58)xs=e−d1ps+2u1s+e−τ2ss+2u2s=e−τ1sg11s+g12u1s+e−τ2sg21s+g22u2ss2+f1s+f2.Such that(59)f1=4,f2=4,g11=1,g12=2,g21=1,g22=2.d1=6.84s,d2=1.87s.


#### 6.2.1. Simulation Results

The simulation results in the case of a noise free or a noisy output are tabulated in Tables5 and 6.Table 5
Simulation results with a noise-free output.
TrueEstimatesf143.9954 ± 0.0046f144.0003 ± 0.0003g1111.0043 ± 0.0043g1221.9747± 0.0253g2110.9838 ± 0.0162g2221.9938 ± 0.0062d16.846.8314 ± 0.0086d21.871.6886 ± 0.0014Table 6
Results simulation withSNR=20dB.
TrueEstimatesf144.0009 ± 0.0009f243.9976 ± 0.0024g1110.9977 ± 0.0023g1222.0007 ± 0.0007g2110.9942 ± 0.0058g2222.0017 ± 0.0017d16.846.8392 ± 0.0008d21.871.8735 ± 0.0035From Tables5 and 6, we can see that the proposed method gives optimal estimates in the noise free and the noisy cases [31, 32].


#### 6.2.2. Computational Time

To further demonstrate the effectiveness of the proposed identification method, a comparison is presented in Table7 with the computational time, respectively, of the GSEPNIV and the proposed method in the case of a noisy output.Table 7
Computational time.
GSEPNIVProposed methodComputational time (sec)4.3350811.713794It can be seen that the proposed prediction error method outperforms the GSEPNIV method in terms of the computational time.


## 6.1. Example 1

The first simulated example is described as(54)x¨t+f1x˙t+f2xt=g11u˙1t−d1+g12u1t−d1+g21u˙2t−d2+g22u2t−d2.In the Laplace domain, the above equation becomes(55)xs=e−d1ss+1u1s+2e−d2ss+2u2s=e−d1ss+2u1s+e−d2s2s+2u2ss2+3s+2.The real parameters are the following:(56)f1=3,f2=2,g11=1,g12=2,g21=2,g22=2.d1=8.83s,d2=2.32s.


### 6.1.1. Simulation Results

We obtain the following results illustrated in Table1.Table 1
Simulation results with noise-free output.
TrueEstimatesf132.9818 ± 0.0182f221.9844 ± 0.0156g1110.9912 ± 0.0088g1221.9876 ± 0.0124g2121.9941 ± 0.0059g2221.9830 ± 0.017d18.838.8301 ± 0.0001d22.322.3255 ± 0.0055With aSNR=20dB, the parameters estimated are shown in Table 2.Table 2
Results simulation withSNR=20dB.
TrueEstimatesf132.9813 ± 0.0187f221.9840 ± 0.016g1110.9931 ± 0.0069g1221.9887 ± 0.0113g2121.9955 ± 0.0045g2221.9810 ± 0.019d18.838.8276 ± 0.0024d22.322.3240 ± 0.0032


### 6.1.2. Model Validation

Validations of the obtained model by using the gradient algorithm are plotted, respectively, in the following Figures2 and 3.Figure 2
Model validation of the proposed scheme with noise-free output.Figure 3
Model validation in the case ofSNR=20dB.


### 6.1.3. Bode Diagrams

The bode diagrams are plotted, respectively, in Figures4 and 5.Figure 4
Bode diagrams of the true system model with the estimated model: noise-free case.Figure 5
Bode diagrams of the true system model with the estimated modelSNR=20dB.The accuracy of the developed prediction error method is gradually improved in the frequency domain.


### 6.1.4. Computational Time

A comparison is done in Table3 with the computational time in the second simulation case. For comparison, the identification algorithm GSEPNIV for the global separable nonlinear instrumental variable given in reference [23] is also adopted to estimate the above stable-type model parameters defined by (53) by taking the same initial parameter setting.Table 3
Computational timeSNR=20dB.
GSEPNIVProposed methodComputational time (sec)5.7307581.756012The computational time is measured using the tic and toc Matlab function.Table3 shows that the GSEPNIV method needs a long computational time. However, it can be seen that the computational time is smaller by using the proposed prediction error method.


### 6.1.5. Noise Effects

Table4 addresses another comparison between the proposed method and the GSEPNIV method developed in [23] in terms of robustness with regard to noise in the case of a noisy output.Table 4
Parameters estimates: simulation 2.
TrueGSEPNIVProposed methodf133.0646 ± 0.06462.9813 ± 0.0187f222.0773 ± 0.07731.9840 ± 0.016g1110.9804 ± 0.01960.9931 ± 0.0069g1222.0759 ± 0.07591.9887 ± 0.0113g2122.0269 ± 0.02691.9955 ± 0.0045g2222.0758 ± 0.07581.9810 ± 0.019d18.848.8671 ± 0.02718.8276 ± 0.0024d22.322.3315 ± 0.01152.3240 ± 0.0032Judging from the simulation results tabulated in Table4 concerning computational error values for each parameter, it can be seen that the proposed method outperforms the GSEPNIV method.


## 6.1.1. Simulation Results

We obtain the following results illustrated in Table1.Table 1
Simulation results with noise-free output.
TrueEstimatesf132.9818 ± 0.0182f221.9844 ± 0.0156g1110.9912 ± 0.0088g1221.9876 ± 0.0124g2121.9941 ± 0.0059g2221.9830 ± 0.017d18.838.8301 ± 0.0001d22.322.3255 ± 0.0055With aSNR=20dB, the parameters estimated are shown in Table 2.Table 2
Results simulation withSNR=20dB.
TrueEstimatesf132.9813 ± 0.0187f221.9840 ± 0.016g1110.9931 ± 0.0069g1221.9887 ± 0.0113g2121.9955 ± 0.0045g2221.9810 ± 0.019d18.838.8276 ± 0.0024d22.322.3240 ± 0.0032


## 6.1.2. Model Validation

Validations of the obtained model by using the gradient algorithm are plotted, respectively, in the following Figures2 and 3.Figure 2
Model validation of the proposed scheme with noise-free output.Figure 3
Model validation in the case ofSNR=20dB.


## 6.1.3. Bode Diagrams

The bode diagrams are plotted, respectively, in Figures4 and 5.Figure 4
Bode diagrams of the true system model with the estimated model: noise-free case.Figure 5
Bode diagrams of the true system model with the estimated modelSNR=20dB.The accuracy of the developed prediction error method is gradually improved in the frequency domain.


## 6.1.4. Computational Time

A comparison is done in Table3 with the computational time in the second simulation case. For comparison, the identification algorithm GSEPNIV for the global separable nonlinear instrumental variable given in reference [23] is also adopted to estimate the above stable-type model parameters defined by (53) by taking the same initial parameter setting.Table 3
Computational timeSNR=20dB.
GSEPNIVProposed methodComputational time (sec)5.7307581.756012The computational time is measured using the tic and toc Matlab function.Table3 shows that the GSEPNIV method needs a long computational time. However, it can be seen that the computational time is smaller by using the proposed prediction error method.


## 6.1.5. Noise Effects

Table4 addresses another comparison between the proposed method and the GSEPNIV method developed in [23] in terms of robustness with regard to noise in the case of a noisy output.Table 4
Parameters estimates: simulation 2.
TrueGSEPNIVProposed methodf133.0646 ± 0.06462.9813 ± 0.0187f222.0773 ± 0.07731.9840 ± 0.016g1110.9804 ± 0.01960.9931 ± 0.0069g1222.0759 ± 0.07591.9887 ± 0.0113g2122.0269 ± 0.02691.9955 ± 0.0045g2222.0758 ± 0.07581.9810 ± 0.019d18.848.8671 ± 0.02718.8276 ± 0.0024d22.322.3315 ± 0.01152.3240 ± 0.0032Judging from the simulation results tabulated in Table4 concerning computational error values for each parameter, it can be seen that the proposed method outperforms the GSEPNIV method.


## 6.2. Example 2

We consider the following process described in Figure6, where u1 is the steam flow, u2 is the water flow u2, and x is the water temperature.Figure 6
The steam-water heat exchanger.The equation of the above process is(57)x¨t+f1x˙t+f2xt=g11u˙1t−d1+g12u1t−d1+g21u˙2t−d2+g22u2t−d2.Then,(58)xs=e−d1ps+2u1s+e−τ2ss+2u2s=e−τ1sg11s+g12u1s+e−τ2sg21s+g22u2ss2+f1s+f2.Such that(59)f1=4,f2=4,g11=1,g12=2,g21=1,g22=2.d1=6.84s,d2=1.87s.


### 6.2.1. Simulation Results

The simulation results in the case of a noise free or a noisy output are tabulated in Tables5 and 6.Table 5
Simulation results with a noise-free output.
TrueEstimatesf143.9954 ± 0.0046f144.0003 ± 0.0003g1111.0043 ± 0.0043g1221.9747± 0.0253g2110.9838 ± 0.0162g2221.9938 ± 0.0062d16.846.8314 ± 0.0086d21.871.6886 ± 0.0014Table 6
Results simulation withSNR=20dB.
TrueEstimatesf144.0009 ± 0.0009f243.9976 ± 0.0024g1110.9977 ± 0.0023g1222.0007 ± 0.0007g2110.9942 ± 0.0058g2222.0017 ± 0.0017d16.846.8392 ± 0.0008d21.871.8735 ± 0.0035From Tables5 and 6, we can see that the proposed method gives optimal estimates in the noise free and the noisy cases [31, 32].


### 6.2.2. Computational Time

To further demonstrate the effectiveness of the proposed identification method, a comparison is presented in Table7 with the computational time, respectively, of the GSEPNIV and the proposed method in the case of a noisy output.Table 7
Computational time.
GSEPNIVProposed methodComputational time (sec)4.3350811.713794It can be seen that the proposed prediction error method outperforms the GSEPNIV method in terms of the computational time.


## 6.2.1. Simulation Results

The simulation results in the case of a noise free or a noisy output are tabulated in Tables5 and 6.Table 5
Simulation results with a noise-free output.
TrueEstimatesf143.9954 ± 0.0046f144.0003 ± 0.0003g1111.0043 ± 0.0043g1221.9747± 0.0253g2110.9838 ± 0.0162g2221.9938 ± 0.0062d16.846.8314 ± 0.0086d21.871.6886 ± 0.0014Table 6
Results simulation withSNR=20dB.
TrueEstimatesf144.0009 ± 0.0009f243.9976 ± 0.0024g1110.9977 ± 0.0023g1222.0007 ± 0.0007g2110.9942 ± 0.0058g2222.0017 ± 0.0017d16.846.8392 ± 0.0008d21.871.8735 ± 0.0035From Tables5 and 6, we can see that the proposed method gives optimal estimates in the noise free and the noisy cases [31, 32].


## 6.2.2. Computational Time

To further demonstrate the effectiveness of the proposed identification method, a comparison is presented in Table7 with the computational time, respectively, of the GSEPNIV and the proposed method in the case of a noisy output.Table 7
Computational time.
GSEPNIVProposed methodComputational time (sec)4.3350811.713794It can be seen that the proposed prediction error method outperforms the GSEPNIV method in terms of the computational time.


## 7. Discussions

From Tables1–7 and Figures 2–5, we can draw the following features:(1)
The proposed method insures a high accuracy in the noise-free case or in the noisy case(2)
According to Figures2–3, we can see clearly that the estimated outputs track fast the true outputs with acceptable speed(3)
The simulation results tabulated in Tables1–7 confirm the robustness of the method proposed in this paper in the case of noise-free output likewise in the case of noisy output(4)
Considering the computational time, the proposed method converges faster than the GSEPNIV method proposed in [23]


## 8. Conclusion

In the present research work, we have addressed the problem of simultaneous identification of the MISO CT system with multiple unknown time delays from sampled data. Indeed, we have developed a new scheme based on a gradient algorithm. The proposed algorithm is then simulated on two examples to confirm its performance.By using the presented scheme, the identification of the considered system operating on the closed loop is under investigation.



---

*Source: 5029272-2022-11-17.xml*