Identification of LPV partial differential equation models

Identification of LPV partial differential equation models Julien Schorsch, Marion Gilson, Vincent Laurain, Hugues Garnier To cite this version: Julien Schorsch, Marion Gilson, Vincent Laurain, Hugues Garnier. Identification of LPV partial differential equation models. 52nd IEEE Conference on Decision and Control, CDC 2013, Dec 2013, Florence, Italy. Paper ThA17.3. ฀hal-00861956฀ HAL Id: hal-00861956 https://hal.archives-ouvertes.fr/hal-00861956 Submitted on 14 Sep 2013 HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers. L’archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d’enseignement et de recherche français ou étrangers, des laboratoires publics ou privés. Identification of LPV partial differential equation models J. Schorsch, M. Gilson, V. Laurain and H. Garnier Abstract— This paper deals with the identification of linear parameter varying (LPV) models described by partial differential equations (PDE). A direct identification of continuous space-time LPV-EDP systems in an input-output setting is investigated in the case of an additive output noise. The continuous space-time LPV-PDE model is firstly proposed to be rewritten as a multi-input single-output linear time-space invariant model and an iterative optimization is then developed to estimate efficiently the model parameters. The performance of the proposed method is then illustrated via a representative simulation example based on an Alsace river-flow measurement. I. I NTRODUCTION Partial differential equations (PDE) describe a wide range of physical behaviors with many applications in various fields from the state-of-health in advanced batteries [1] to larger scale problems such as simulation of the Arctic sea ice temperature [2]. The main characteristics of PDE with respect to the usually considered ordinary differential equations (ODE) is the introduction of partial derivatives along different dimensions, mostly depicting the time and space components of the considered system. PDE models are most of the time derived from physical laws, and therefore expressed in continuous time-space (CTS). To solve the identification problem of PDE systems, two main approaches can be distinguished: • Indirect approaches which includes two steps. The first involves the identification of a discrete time-space (DTS) model. In the second step, the identified DTS model is converted into CTS (see e.g. [13]). • Direct approaches: which aim at estimating the coefficients of the CTS model directly from sampled measured signals (see e.g. [15], [14], [16]). In the usual ODE linear system case, one of the main limitation of indirect approaches is the impossibility to deal with non-regularly sampled data. In the presented PDE case, while it is fair to assume a regular time-sampling, the spatial-sampling is directly bounded to a physical sensor implantation and hence limited by cost and physical possible positioning. Consequently, in practice, indirect identification of PDE systems might often not be applicable. In order to cope with these limitations, direct approaches such as proposed in [15] aim at directly identifying the physical CTS parameters, with the inherent difficulty to approximate J. Schorsch, M. Gilson, V. Laurain and H. Garnier are with University of Lorraine, CRAN, UMR 7039, 2 rue Jean Lamour, F-54519 Vandœuvre-les-Nancy, France, CNRS, CRAN, UMR 7039, France julien.schorsch, marion.gilson, vincent.laurain, hugues.garnier@univ-lorraine.fr partial derivatives from sampled noisy signals. This difficulty is commonly solved by a low pass filtering of the data. Nevertheless, it was only recently that the noise modeling problem with respect to these filters was addressed for direct identification of PDE systems (see e.g. [7], [19]). Another statement is that the PDE considered for system identification are mostly linear time and space invariant (LTSI). Nevertheless, in practice most of the physical behaviors have a nonlinear or time-space varying nature. Consider for example a pollutant transport along a river. From a conceptual reasoning, it clearly appears that the transport speed directly depends on the river flow or the river profile for example. Consequently, this paper aims at the identification of linear parameter varying (LPV) PDE systems. Identification of LPV systems has gained some strong interest in the past few years, mostly in the discretetime (DT) case [10], [11]. Identification of DTS-LPV-PDE has been dealt with in [3], using an instrumental variable (IV) method. The use of LPV models has also been applied to the identification of LTSI-PDE system by approximating the PDE at each working point with a DT-LPV model such as in [4], [5]. Regarding direct identification of CT-LPV models, only a few methods are available in the ODE case [6], and to the best of authors’ knowledge, there is no specific work in the PDE case. Consequently, this paper aims at identifying CTSLPV-PDE systems from sampled data under some realistic noise conditions. As pointed out in [6], the identification of CT parameters with physical meaning in LPV systems requires the use of a direct identification scheme which offers challenging theoretical issues. In the DT ODE case, a convolution type representation of the system dynamics, i.e. an LPV Impulse Response Representation (IRR) [10] is used in order to define predictors with respect to some model class. In the CTS case however, no IRR has been developed yet and thus the same concept cannot be used to define the predictors such as usually done [12]. Yet, a system reformulation was proposed in [6] in order to apply the linear time invariant (LTI) system theory to the LPV identification problem. This reformulation is here extended to the PDE case, allowing the use of a least squares method which provides efficient estimates for output-error type of noise in the prediction error minimization (PEM) sense. The paper is organized as follows: in Section 2, the general class of CTS-LPV-PDE systems in an input-output (IO) representation intended for identification is introduced, along with the main issues related to direct CTS identification. In section 3, the reformulation of the system which enables the predictor quality definition is provided. Finally, an iterative least-squares (LS) algorithm dedicated to this model identification is detailed in section 4 and it is thoroughly analyzed on a relevant example in Section 5. II. P ROBLEM FORMULATION A. System description Consider the data generating CTS-LPV system described by the following spatio-temporal PDE: S   Ao (p, ρ)χ̊(t, x) = Bo (p, ρ)u(t, x) (1)  B. Model considered It can be seen that the analysed process is fully characterized by the knowledge of functions {åit ,ix } and {b̊jt ,jx }. Yet a main difficulty in the estimation of LPV systems is that in practical cases, these functions are a priori unknown nonlinear functions. Even if some recent work gives some preliminary results oriented towards a model structure selection, the usual solution to overcome this problem is to parameterize these functions using a sum of a priori known basis functions [6], [10]. Consequently, this paper focuses on the usual assumption used in the literature where the process model denoted by Gθ is defined in the following form: y(t, x) = χ̊(t, x) + vo (t, x) Gθ : This equation is subject to initial conditions and/or boundary conditions. u(t, x) denotes the input variable observed at time t and at space point x, χ̊(t, x) the noise-free output and vo (t, x) the additive output noise. The parameter varying nature of the system behavior is caracterized by ρ, the socalled scheduling variable with ρ = ρ(t, x) denoting the dependency of the system on external deterministic signal both in time and space. It can be noticed that as a first step towards more general noise assumptions, throughout this paper, the noise is supposed to be independent on ρ. In this PDE context, p = (pt , px ) denotes the partial differential operators: pit ,ix u(t, x) = pitt pixx u(t, x) = ∂ it +ix u(t, x) ∂tit ∂xix (2) ∂ ∂x . ∂ ∂t and px = Ao (p, ρ) and Bo (p, ρ) are where pt = bivariate polynomials in p of degrees nt , nx and mt , mx respectively: (B(p, ρ, θ), A(p, ρ, θ)) = (Aθ , Bθ ) (4) where the ρ-dependent polynomials A and B are parameterized as  nt ,nx  (it ,ix )6=(nt ,nx )  X    ait ,ix (ρ)pit ,ix A(p, ρ, θ) = pnt ,nx +     it ,ix =0  Aθ : (5) nα  X  a 0  alit ,ix fl (ρ)  i ,i (ρ) = ait ,ix +   t x  l=1   it = 0, · · · , nt ix = 0, · · · , nx Bθ :  mX t ,mx    bjt ,jx (ρ)pjt ,jx B(p, ρ, θ) =     jt ,jx =0   mβ X   0  + bljt ,jx gl (ρ) (ρ) = b b  jt ,jx jt ,jx    l=1   jt = 0, · · · , mt jx = 0, · · · , mx (6) m     nt ,nx  +   Ao (p, ρ) = p nt ,nx (it ,ix )6=(nt ,nx ) X åit ,ix (ρ)pit ,ix (3a) it ,ix =0 mX  t ,mx     b̊jt ,jx (ρ)pjt ,jx  Bo (p, ρ) = (3b) jt ,jx =0 where åit ,ix (ρ) and b̊jt ,jx (ρ) are analytic functions with dependence on ρ at time t and at space x, which is called static dependence, note that ρ = ρ(t, x). Like in any identification problem, it is assumed that a data sequence DN,L {y(tk , xℓ ), u(tk , xℓ ), ρ(tk , xℓ )}N,L k,ℓ=1 generated by S is available. Let h be the spatial sampling distance between two observation points, and Ts the temporal sampling interval. L denotes the number of spatial observation locations and N the number of temporal sampling instants. The data-generating system in Eq (1) is then identified by defining a model set and some information criterion characterizing the quality of a given model based on DN,L . As a next step, a model class suitable for the identification of such systems is given. β α In this parametrization, {fl }nl=1 and {gl }l=1 are meromorphic functions of ρ, with static dependence. They can be chosen for examples as linearly independent functions. The associated model parameters θ are stacked columnwise: T θ = [a0,0 · · · ant ,nx −1 b0,0 · · · bmt ,mx ] ∈ Rnθ (7) where   (8) ait ,ix = a0it ,ix a1it ,ix · · · anitα,ix ∈ Rnα +1  mβ  ∈ Rmβ +1 bjt ,jx = b0jt ,jx b1jt ,jx · · · bjt ,j x (9) and nθ = ((nt + 1)(nx + 1) − 1)(nα + 1) + (mt + 1)(mx + 1)(mβ + 1). Let introduce also G = {Gθ | θ ∈ Rnθ }, as a collection of all process models in the form of (4). Furthermore, it has to be noticed that the system presented in Eq (1) displays a CTS representation of the noise. Nevertheless, in order to avoid the rather difficult mathematical problem of treating sampled CTS random process in terms of a filtered piece-wise constant CTS noise source (see [9], [8]), it is here proposed to consider only the sampling instances of the noise and that consequently, the model denoted as Mθ is then defined as the following hybrid structure [18]:  A(p, ρ, θ)χ(t, x) = B(p, ρ, θ)u(t, x)  Mθ : (10)  y(tk , xℓ ) = χ(tk , xℓ ) + e(tk , xℓ ) Furthermore, for clarity’s sake, it is considered that the additive noise e(xℓ , tk ) is a zero mean, random spatial array, with no temporal or spatial correlation and that vo (tk , xℓ ) also fulfills these assumptions. In order to handle the more complicated case of a filtered DT noise, some hints can be found in [6]. Even though the system behavior is completely characterized by the LTSI vector parameter θ in the model structure defined in (4), the quality assessment of a given estimator can not be straight-forwardly defined in the CTS-LPV setting. A solution to this problem is developed in the next section where a system reformulation based on [6] is proposed, allowing the definition of a prediction error minimization framework in the CTS-LVP-PDE case. III. R EFORMULATION now expressed as  (i ,i )  χ̊l t x (t, x) = fl (ρ)χ̊(it ,ix ) (t, x) (13a)     {it , ix , l} ∈ {0 · · · nt , 0 · · · nx , 1 · · · nα }  u(jt ,jx ) (t, x) = gl (ρ)u(jt ,jx ) (t, x) (13b)  l    {jt , jx , l} ∈ {0 · · · mt , 0 · · · mx , 0 · · · mβ } Therefore, the process part of the LPV model is rewritten as a MISO system with ((nt + 1)(nx + 1) − 1)nα + ,nx ,nα (mt + 1)(mx + 1)(mβ + 1) inputs {χlit ,ix }intt=0,i and x =0,l=1 mt ,mx ,mβ l {ujt ,jx }jt =0,jx =0,l=0 . By using (11), (1) can be rewritten in terms of the sampled output signal y(tk , xℓ ) as  nt ,nx (it ,ix )6=(nt ,nx ) nα  X X åli ,i (i ,i )   t x   χ̊l t x (tk , xℓ ) y(tk , xℓ ) = − F (p) o it ,ix =0 l=1 (14) mβ l mt ,mx X   X b̊jt ,jx (jt ,jx )   u (tk , xℓ ) + eo (tk , xℓ ).  + Fo (p) l j ,j =0 t x l=0 At this point, it must be noticed that despite the apparent writing complexity, equations (14) and (1) are equivalent. Therefore the model to be identified is still given by equation (10). The gain here is the LTSI formulation of the LPVPDE system which consequently enables the use of a widely spread theory on linear systems for the expression of the optimization in the prediction error sense. The obvious drawback is the lack of a priori knowledge on noise-free outputs χ̊it ,ix . A solution to cope with this issue is proposed in Section IV. However, by momentarily assuming the knowledge of χ̊it ,ix , the prediction error can be expressed as A PEM cannot be directly applied to identify a model described by (10) because of the non-commutativity of the LPV polynomials. This problem is well-known when dealing with LPV models. In [6], the solution proposed to identify a LPV-ODE model is therefore to rewrite the LPV system as a multi-input single-output (MISO) LTI model for which the PEM can be applied. εθ (tk , xℓ ) = y(tk , xℓ ) − ŷθ (tk , xℓ ) (15) The same idea is used here for LPV-EDP model. Indeed, if the system belongs to the model set defined in (10), the where ŷθ (tk , xℓ ) is the one-step-ahead predictor [12]: LPV-PDE model can be expressed as a CTS-MISO model nt ,nx (it ,ix )6=(nt ,nx ) nα X X ali ,i by rewriting the data-generating system (1) as the following (i ,i ) t x χ̊l t x (tk , xℓ ) ŷθ (tk , xℓ ) = − LTSI representation F (p, θ) it ,ix =0 l=1  (16) nt ,nx mβ mX t ,mx X l  (it ,ix )6=(nt ,nx ) b  j ,j (j ,j ) X t x t x  u (tk , xℓ ) . +  (nt ,nx )  (t, x) + å0it ,ix χ̊(it ,ix ) (t, x)  F (p, θ) l χ̊ j ,j =0 l=0  t x  it ,ix =0  | {z }  Based on these considerations, the following section   Fo (p)χ̊(t,x)   presents a LS-based identification method which provides nt ,nx   (it ,ix )6=(nt ,nx ) nα X X efficient estimates i.e: + ålit ,ix fl (ρ)χ̊(it ,ix ) (t, x) (11)  N,L  | {z }  1 X 2  it ,ix =0 l=1  (it ,ix ) θ̂ = arg min εθ (tk , xℓ ), (17)  χ̊ (t,x)  l  θ∈Rnθ N L mβ  mX t ,mx X k,ℓ=1     = b̊ljt ,jx gl (ρ)u(jt ,jx )(t,x)  for data generating systems S written as in (1) and such that  {z } |   jt ,jx =0 l=0  (j ,j ) S ∈ Mθ is informative with respect to Mθ . u t x (t,x) l where g0 (ρ) = 1 and the superscript (jt , jx ) for a signal, like u(jt ,jx ) , denotes the jtth partial derivative in time and the jxth partial derivative in space and y (nt ,nx ) (tk , xℓ ) = ϕT (tk , xℓ )θo + Fo (p)eo (tk , xℓ ) nt ,nx (it ,ix )6=(nt ,nx ) Fo (p) = pnt ,nx + X IV. I TERATIVE LEAST SQUARES FOR LPV-EDP SYSTEMS Using the LTSI system reformulation (14), y(tk , xℓ ) can be rewritten in the following regression form: å0it ,ix pit ,ix (12) it ,ix =0 Note that in this way, the main advantage is to transpose the time-space variation of the coefficients onto the signals (18) where,  ϕ(tk , xℓ ) = −y(tk , xℓ ) · · · − y (nt ,nx −1) (tk , xℓ ) (0,0) (n ,n −1) −χ̊1 (tk , xℓ ) · · · − χ̊nαt x (tk , xℓ ) (0,0) (mt ,mx ) T u0 (tk , xℓ ) · · · umβ (tk , xℓ ) ] (19) h iT mβ 0 α θo = å00,0 · · ·å0nt ,nx −1 å10,0 · · ·ån (20) nt ,nx −1 b̊0,0 · · · b̊mt ,mx In the given context, contrary to the LPV formulation where the LPV dependency appears into the signal definitions, the filter Fo (p) commutes and allows the rewritten of (18) as follows (nt ,nx ) yf (tk , xℓ ) = ϕT f (tk , xℓ )θo + eo (tk , xℓ ) (i ,i ) (21) (j ,j ) where yf (tk , xℓ ), χ̊l,ft x (tk , xℓ ) and ul,ft x (tk , xℓ ) represent the outputs of a prefiltering operation. 1 y(tk , xℓ ) (22) yf (tk , xℓ ) = Fo (p) and h (n ,n −1) ϕf (tk , xℓ ) = −yf (tk , xℓ ) · · · − yf t x (tk , xℓ ) (0,0) (n ,n −1) −χ̊1,f (tk , xℓ ) · · · − χ̊nαt,f x (tk , xℓ ) (0,0) (m ,m ) T u0,f (tk , xℓ ) · · · umβt,f x (tk , xℓ ) ] (23) It has to be noticed that the interest of the filter Fo (p) is twofold. Firstly it is used to filter the input/output signals used for the identification to come up with the regression equation (21) written in terms of white noise eo (tk , xℓ ). In other words, Fo (p) is used to filter the data such that a simple LS algorithm applied to these pre-filtered data leads to the statistically optimal estimate under the condition S ∈ Mθ . Secondly, another important particularity of this method is to provide directly a CTS model, therefore the regressor ϕf (tk , xℓ ) (23) is built up form the time-space derivatives of the data (and not from the delayed version of these data as in the DTS case). The main problem to be solved is that these derivatives are unknown and has therefore to be estimated. This well-known problem is usually solved by using low pass filtering on the input and output signals, as e.g. statevariable-filters. However these filters present the drawback of requiring the critical choice of a design variable. On the contrary, the particular advantage of the presented method is to be free of this drawback since the filter Fo (p) is not only used for the prediction error minimization but also for the approximation of the unknown derivatives. As a result, the filtered regressor ϕf (tk , xℓ ) in (21) is directly numerically approximated instead of the regressor ϕ(tk , xℓ ) in (18) which makes use of unknown time-space derivatives and still requires a pre-filtered operation to access the space-time derivatives, with a filter which has to be chosen by the user. Consequently, under the assumption of a priori known Fo (p) and χ̊(it ,ix ) , the optimal estimator in the prediction error minimization sense is the LS estimator  −1 N,L X  · θ̂LS =  ϕf (tk , xℓ )ϕT f (tk , xℓ ) k,ℓ=1   N,L X k,ℓ=1  (24) (n ,n ) ϕf (tk , xℓ )yf t x (tk , xℓ ) (it ,ix ) Naturally, the filter Fo (p) and the signals χ̊ are unknown and therefore and iterative procedure is proposed in order to cope with this issue, similarly as it is achieved in the approaches developed for the ODE case in [21], [22]. From Eq. (24) two comments may be formulated. Indeed the identification problem defined by Eq. (21) has been chosen to be solved by an iterative LS method. Alternatives could have also been chosen, as e.g. more advanced PEM method [12]. However, these latter need a robust initialization which seems not easy to handle in this CTS-PDE-LPV model identification case. Moreover, the iterative solution appears as a natural choice since the proposed method involves an iterative solution (or relaxation) algorithm in which, at each iteration, an auxiliary model is used to generate a part of the regressor ϕf (tk , xℓ ) as well as the associated pre-filters. As a result the regressor used for the identification is not only based on reconstructed signals but on real input, measured output, and reconstructed output. LPV-LS-PDE algorithm: The proposed LPV-LS-PDE model identification algorithm can be summarized as shown below. The algorithm consists of three first steps for the initialization and five iterative steps used to estimate the filter, the derivatives and finally the model parameters. • • Step 1 Initialization of θ0 for l = 0 only: iT h θ̂0 = â00,0 · · · â0nt ,nx −1 0 · · · 0b̂00,0 · · · b̂0mt ,mx 0 · · · 0 using for instance the LS-PDE method proposed in [19]. In this first step the system is considered as a LTSI system. Step 2 Compute a first estimated continuous-filter on the basis of the estimates obtained in Step 1 1 Q(p, θ̂0 ) = F (p, θ̂0 ) where F (p, θ̂0 ) is given as (12). Use the filter in order to generate the estimates of the derivatives (j ,jx ) mt ,mx ,mβ }jt =0,jx =0,l=0 {ul,ft (it ,ix ) nt ,nx }it =0,ix =0 {yf (jt ,jx ) mt ,mx ,mβ }jt =0,jx =0,l=0 = Q(p, θ̂0 ){ul ,nx = Q(p, θ̂0 ){y (it ,ix ) }intt=0,i x =0 and the regressor ϕf (tk , xℓ ) built up from the filtered noisy output only and the input for this initialization step h (n ,n −1) ϕf (tk , xℓ ) = −yf (tk , xℓ ) · · · − yf t x (tk , xℓ ) (0,0) (n ,n −1) −y1,f (tk , xℓ ) · · · − ynαt,f x (tk , xℓ ) (0,0) (m ,m ) T u0,f (tk , xℓ ) · · · umβt,f x (tk , xℓ ) ] • step 3 Compute an estimate θ̂1 : θ̂1 =   −1 N,L X k,ℓ=1  ϕf (tk , xℓ )ϕT f (tk , xℓ )   · N,L X k,ℓ=1 (nt ,nx ) ϕf (tk , xℓ )yf  (tk , xℓ ) Consequently, the initial estimate proposed for the LPV-LS-EDP algorithm is an LTSI-LS-EDP estimate of Mθ . Set τ = 1. [End of initialization part.] analysis, the boundary and initial conditions are assumed to be zeros. V. E XAMPLE • Step 4 At any iteration τ +1, compute an estimate of the ˆ k , xℓ ) via numerical approximation noise-free output χ̊(t of ˆ x) = B(p, ρ, θ̂τ )u(t, x) A(p, ρ, θ̂τ )χ̊(t, based on the estimates θ̂τ at the previous iteration. Step 5 Compute the estimated continuous filter Q(p, θ̂τ ) = 1 F (p, θ̂τ ) ˆ k , xℓ ) in order to generate and use the filter as well as χ̊(t the estimates of the derivatives which are needed: (it ,ix ) nt ,nx ,nα }it =0,ix =0,l=1 ˆ (it ,ix ) }nt ,nx ,nα = Q(p, θ̂τ ){χ̊ l it =0,ix =0,l=1 (j ,jx ) mt ,mx ,mβ }jt =0,jx =0,l=0 = Q(p, θ̂τ ){ul ˆ {χ̊ l,f {ul,ft (it ,ix ) nt ,nx }it =0,ix =0 {yf • (jt ,jx ) mt ,mx ,mβ }jt =0,jx =0,l=0 ,nx = Q(p, θ̂τ ){y (it ,ix ) }intt=0,i x =0 Step 6 Build the filtered estimated regressor ϕ̂f (tk , xℓ ) as: h (n ,n −1) ϕ̂f (tk , xℓ ) = −yf (tk , xℓ ) · · · − yf t x (tk , xℓ ) ˆ (0,0) (tk , xℓ ) · · · − χ̊ ˆ (nt ,nx −1) (tk , xℓ ) −χ̊ 1,f nα ,f (0,0) (m ,m ) T u0,f (tk , xℓ ) · · · umβt,f x (tk , xℓ ) ] • ϕ̂f (tk , xℓ ) vector is composed of the filtered noisy output, the filtered noise-free estimated output and the input Step 7 Compute the LPV-LS-PDE estimate θ̂τ +1 at the (τ + 1)th iteration: θ̂ τ +1 =   −1 N,L X k,ℓ=1  ϕ̂f (tk , xℓ )ϕ̂T f (tk , xℓ )   • N,L X k,ℓ=1 ·  (n ,n ) ϕ̂f (tk , xℓ )yf t x (tk , xℓ ) Step 8 If θ̂τ +1 has converged or the maximum number of iterations is reached, then stop, else increase τ by 1 and go to Step 4. Comments on the implementation issue: In this paper, we considered the practical feasible situation such that only sampled measurements of the spacetime continuous signals (y, u, ρ) are available. In order to apply a continuous filter on sampled data one can either interpolate the samples to obtain a continuous signal and apply the continuous filter on this reconstructed signal or use a numerical approximation. This is a common problem for simulation of continuous systems. For simulation purposes, discrete approximation of the system can efficiently be dealt with by using available powerful numerical algorithms. Note in the following example, an Euler method is used to simulate the true system and the filtering operations. In the above The performance of the described algorithm is presented on a representative simulation example. As an example, the advection-diffusion equation (ADE) is used. This equation is often used in the water resource quality analysis to describe the transport and dispersion of a solute (pesticide, pollutant, ...) into a river channel (see e.g. [20]). The ADE can be written as following:  ∂ 2 χ̊(t, x) ∂χ̊(t, x)   = å02 (ρ(t))    ∂x2  ∂t ∂χ̊(t, x) (25) −å01 (ρ(t)) + b̊00 (ρ(t))u(t, x)  ∂x   χ̊(t = 0, x) = 0    y(tk , xℓ ) = χ̊(tk , xℓ ) + eo (tk , xℓ ) where the noise-free signal χ̊(xℓ , tk ) (Fig. 1) represents the solute propagation in the river and is obtained by discretization of the partial differential equation by a finite difference method (see e.g. [17]). the temporal sampling is Ts = 1day and N = 3650 temporal points are measured. The spatial sampling is h = 20m for L = 50 space points. These 50 points can represent 50 sensors distributed in the river. u(t, x) describes the source of the pollutant, and is represented by a PRBS signal at the space point x = 20m. In realistic conditions of simulation the a02 (ρ(t)) diffusion 0 18 100 16 200 14 300 distance in meter • 12 400 10 500 8 600 700 6 800 4 900 2 1000 0 0 500 1000 1500 2000 2500 3000 3500 Time in day Fig. 1. Noise-free output χ̊(t, x) and the a01 (ρ(t)) advection coefficients are dependent of the river flow Q(t) (Fig. 2). The scheduling variable ρ is then taken as ρ(t) = Q(t) and the coefficients a02 (ρ(t)) , a01 (ρ(t)) and b00 (ρ(t)) are given by  (26a)   å02 (ρ(t)) = 2 + 0.9ρ(t) å01 (ρ(t)) = 0.5 + 0.2ρ(t) (26b)   b̊00 (ρ(t)) = 1 (26c) In this example, realistic data are used since the river flow considered here comes from measurements of a river located in Alsace (France). 10 years of measurements are used in this example. For model identification purposes, the output signal is supposed to be corrupted by a two dimensional, zero-mean R EFERENCES 120 100 80 60 40 20 0 0 500 1000 1500 2000 2500 3000 3500 Time in day Fig. 2. Flow river Q(t) in m3 /s and normally distributed, discrete-time noise signal, with a Signal-to-Noise Ratio (SNR) defined as   Pχ̊ (27) SN R = 10log10 Pe where Pe and Pχ̊ are the average power of the noise and deterministic output.To provide representative results, a Monte Carlo simulation analysis for 100 runs is used to illustrate the performance of the approach. In Table I and II, it can be noted that the estimates of the proposed method are unbiased with a relatively small standard deviation when SN R = 30dB. Moreover, in the case where SN R = 15dB, the estimates are again unbiased with a higher standard deviation. TABLE I MC SIMULATION RESULTS FOR 100 RUNS ( DIFFUSION COEFFICIENTS ) True value SNR mean std â002 2 30dB 2.0218 0.1071 â102 0.9 15dB 2.0149 0.43073 30dB 0.8922 0.0446 15dB 0.8695 0.34926 TABLE II MC SIMULATION RESULTS FOR 100 RUNS ( ADVECTION COEFFICIENTS ) True value SNR mean std â001 0.5 30dB 0.4993 0.0037 â101 0.2 15dB 0.4998 0.0165 30dB 0.2002 0.0014 b̂000 1 15dB 0.2002 0.0059 30dB 1.0003 0.0016 15dB 0.9981 0.0195 VI. C ONCLUSION A method has been proposed in this paper to identify a continuous space-time LPV partial differential equation model. The approach is based on a MISO linear time/space invariant reformulation of the data generating system along with an iterative procedure develop to estimate efficiently the model parameters into the error-prediction minimization framework. A realistic example has been used to illustrate the efficiency of the proposed method. Future research investigates the use of an instrumental variable method to handle the situation where the additive noise is colored. [1] S. J. Moura, N. A. Chaturvedi, and M. Krstić, ”PDE estimation techniques for advanced battery management systems - Part I: SOC estimation”, Proceedings of the American Control Conference, Montreal (QC), pp: 559 - 565, 2012. [2] H. Fang, J. Wang, E. Feng, and Z. Li, ”Parameter identification and application of a distributed parameter coupled system with a movable inner boundary”, Computers and Mathematics with Applications, vol. 62, pp: 4015-4020, 2011. [3] M. Ali, A. Ali, H. Abbas, and H. Werner, ”Identification of BoxJenkins models for parameter-varying spatially interconnected systems”, American Control Conference, San Francisco (USA), pp: 145150, 2011. [4] G. Belforte, F. Dabbene, and P. Gay, ”LPV approximation of distributed parameter systems in environmental modelling”, Environmental Modelling & Software, vol. 20, pp: 1063-1070, 2005. [5] V. Puig, J. Quevedo, T. Escobet, P. Charbonnaud, and E. Duviella, ”Identification and control of an open-flow canal using LPV models”, Proceedings of the 44th IEEE Conference on Decision and Control, and the European Control Conference, Seville (Spain), pp: 1893-1898, 2005. [6] V. Laurain, R. Tóth, M. Gilson, and H. Garnier, ” Direct identification of continuous-time linear parameter-varying input/output models”, IET Control Theory & Applications, vol. 5(7), pp: 878-888, 2011. [7] J. Schorsch, H. Garnier, and M. Gilson, ”Instrumental variable methods for identifying partial differential equation models of distributed parameter systems”, In IFAC Symposium on System Identification, pp: 840–845, Brussels (Belgium), 2012. [8] R. Pintelon, J. Schoukens, and Y. Rolain, ”Box-Jenkins continuoustime modeling”, Automatica, vol. 36(7), pp: 983-991, 2000. [9] R. Johansson, ’Identification of continuous-time models”, IEEE Transactions on Signal Processing, vol. 42(4), pp: 3887-3897, 1994. [10] R. Tóth. Modeling and Identification of Linear Parameter-Varying Systems. Lecture Notes in Control and Information Sciences, Vol. 403. Springer-Germany, 2010. [11] V. Laurain and M. Gilson and R. Tóth and H. Garnier, ”Refined Instrumental Variable Methods for Identification of LPV Box-Jenkins Models”, Automatica, vol. 46(6), pp: 959-967, 2010. [12] L. Ljung. System Identification, theory for the user. Prentice Hall, 1999. [13] M. Ali and S. S. Chughtai, and H. Werner, ”Consistent identification of two-dimensional systems”, American Control Conference, Baltimore USA, pp: 3464-3469, 2010. [14] M. S. Sadabadi, M. Shafiee, and M. Karrari, ”System identification of two dimensional continuous-time systems using wavelets as modulating functions”, ISA Transactions, vol. 47, pp: 256-266, 2008. [15] C. Chochol, S. Chesne, and D. Remond, ”An original differentiation tool for identification on continuous structures”, Journal of Sound and Vibration, 2013. [16] S. Sagara, Z. J. Yang, and K. Wada, ”Parameter identification of distributed-parameter systems in the presence of measurement noise”, Int. Journal of Systems Science, vol. 22(8), pp: 13911401, 1991. [17] K. W. Morton and D. F. Mayers, Numerical Solution of Partial Differential Equations, Cambridge University Press, 2005. [18] H. Garnier and L. Wang (Eds.), editor, Identification of continuoustime models from sampled data, Springer-Verlag, London, 2008. [19] J. Schorsch, H. Garnier, M. Gilson, and P. C. Young, ”Instrumental variable methods for identifying partial differential equation models”, International Journal of Control, in press, 2013. [20] W. Tych and P.C. Young, ”A matlab software framework for dynamic model emulation”, Environmental Modelling and Software, vol. 34, pp: 19-29, 2012. [21] P. Lohnberg and G. H. J. Wisselink, ”Iterative Least-Squares Parameter Estimation for ARMA Pulse Response and Output Disturbance”, IEEE Transactions on Automatique Control, 27(6), pp: 1252-1255, 1982. [22] F. Dinga, Y. Shi, and T. Chen, ”Auxiliary model-based least-squares identification methods for Hammerstein output-error systems”, Systems & Control Letters, vol. 56, pp: 373 380, 2007.