C Casenave
email: casenave@laas.fr

Identification of Time-non local Models under Diffusive Representation

Keywords: Least squares identification, Sampled data, Operators, Convolution integral, Continuous time systems, Dynamic systems, Distributed-parameter systems

come

INTRODUCTION

Numerous models of physics involve dynamic convolution operators of the form:

u → H(∂ t )u, [H(∂ t )u](t) = ∫ t 0 h(t -s) u(s) ds, ( 1
)
where h is the impulse response of H(∂ t ). Identifying such an operator means identifying one of its characteristic quantities. As the operator to be identified is linear, a convenient and rather general approach consists in working in the frequency domain where any causal operator can be well-defined by its symbol H(iω), that is the Fourier transform of its impulse response. Then, the problem of identifying H(iω) can be classically solved from physical measurements by means of Fourier techniques. Note however that purely frequency identification presents some well-known shortcomings. In particular, the so-identified symbol H(iω) is in general ill adapted to the construction of efficient time-realizations. This is partly due to excessive numerical cost of quadrature approximations resulting from the intrinsic convolution nature of the associated operator, sometimes with long memory [START_REF] Rumeau | Behavior modeling of a cacu 3 ti 4 o 12 ceramic for capacitor applications[END_REF]) or even delay-like behaviors [START_REF] Montseny | Diffusive representation for operators involving delays[END_REF]).

Another shortcoming is that frequency methods are incompatible with real-time identification (and so with pursuit when the symbol has the ability to evolve slowly). But above all, the number of unknown parameters is excessive, what makes the problem excessively sensitive to measurement noises.

In opposite, time domain identification techniques do not present such drawbacks. However, their scope is in general not so wide. See for example [START_REF] Monin | ARMA lattice identification: a new hereditary algorithm[END_REF]) for an interesting optimal method based on ARMA lattices.

By means of diffusive representation, it is possible to use the notion of symbol in the time-domain to identify a small number of unknown parameters, sufficient to get a good accuracy. Indeed, the diffusive representation theory is devoted to state realizations of causal integral operators which can be characterized by a suitable distribution called the γ-symbol. Instead of identifying the symbol H(iω) of the operator, we then identify its γ-symbol µ, which is linearly involved.

Identifying the γ-symbol of a convolution operator presents numerous advantages:

• an input-output stable differential formulation;

• recursive identification algorithms, compatible with real-time identification or even pursuit, can be build; • as for frequency methods, the rational and non rational operators are considered without distinction and can be identified by the same identification process.

The paper is organized as follows. In section 2, we present a simplified version of the diffusive representation. In section 3, we describe the method of identification of a dynamic operator by means of its γ-symbol. Then in section 4, we present two methods of identification of Volterra models. Finally, in section 5, we implement the methods described on a concrete example.

DIFFUSIVE FORMULATION OF CAUSAL CONVOLUTION OPERATORS

A complete statement of diffusive representation will be found in [START_REF] Montseny | Représentation diffusive[END_REF]). Various applications and questions relating to this approach will be found for example in [START_REF] Audounet | A threshold phenomenon in the propagation of a point source initiated flame[END_REF], [START_REF] Carmona | Fractional brownian motion and the markov property[END_REF], Degerli et al. [Ap. 1999], [START_REF] Garcia | Identification of the dynamics of a lead acid battery by a diffusive model[END_REF], [START_REF] Lenczner | Diffusive realization of operator solutions of certain operational partial differential equations[END_REF], [START_REF] Levadoux | Diffusive formulation of the impedance operator on circular boundary for 2D wave equation[END_REF], [START_REF] Montseny | Diffusive representation for operators involving delays[END_REF], [START_REF] Mouyon | Identification of a 2d turbulent wind spectrum[END_REF], [START_REF] Rumeau | Behavior modeling of a cacu 3 ti 4 o 12 ceramic for capacitor applications[END_REF]).

Mathematical framework

We consider a causal convolution operator defined, on any continuous function u : R + → R, by:

u → ( t → ∫ t 0 k(t -s) u(s) ds ) .
(2)

We denote K the Laplace transform of k and K(∂ t ) the convolution operator defined by (2).

Let u t (s) = 1 ]-∞,t] (s) u(s) be the restriction of u to its past and u t (s) = u t (t -s) the so-called "history" of u.

From causality of K(∂ t ), we deduce: [K(∂ t )(u -u t )](t) = 0 for all t;

(3) then, we have for any continuous function u:

[K(∂ t ) u](t) = [ L -1 (K Lu) ] (t) = [ L -1 ( K Lu t )] (t). ( 4
)
We define Ψ u (t, p) := e p t (Lu t ) (p) = (Lu t ) (-p); by computing ∂ t Lu t , Laplace inversion and use of (4): Lemma 1. Ψ u is solution of the differential equation: 5) and for any b 0:

∂ t Ψ(t, p) = p Ψ(t, p) + u, t > 0, Ψ(0, p) = 0, (
[K(∂ t ) u] (t) = 1 2iπ ∫ b+i∞ b-i∞ K(p) Ψ u (t, p) dp. ( 6
)
Let γ be a closed 1 simple arc in C -; we denote Ω + γ the exterior domain defined by γ, and Ω - γ the complementary of Ω + γ . Assume there exists

α γ ∈] π 2 , π[ and a ∈ R such that: e i[-αγ , αγ ] R + + a ⊂ Ω + γ .
(7) By use of standard techniques (Cauchy theorem, Jordan lemma), it can be shown:

Lemma 2. For γ such that K is holomorphic in Ω + γ , if K(p) → 0 when p → ∞ in Ω + γ , then: [K(∂ t ) u] (t) = 1 2iπ ∫ γ K(p) Ψ u (t, p) dp, ( 8
)
where γ is any closed simple arc in Ω + γ such that γ ⊂ Ω - γ .

We now suppose that γ, γ are defined by functions of the Sobolev space 2 W 1,∞ loc (R; C), also denoted γ, γ and such that γ(0) = 0. We use the convenient notation ⟨µ, ψ⟩ = ∫ µ ψ dξ; in particular, when µ is atomic that is

µ = ∑ k a k δ ξ k , we have: ⟨µ, ψ⟩ = ∑ k a k ψ(ξ k ).
Under hypothesis of lemma 2, we have [START_REF] Montseny | Représentation diffusive[END_REF]): Theorem 3. If the possible singularities of K on γ are simple poles or branching points such that |K • γ| is locally integrable in their neighborhood, then:

1. with μ = γ′ 2iπ K • γ and ψ(t, .) = Ψ u (t, .) • γ: [K(∂ t ) u] (t) = ⟨ μ, ψ(t, .) ⟩ ; (9) 2. with 3 γn → γ in W 1,∞ loc and µ = γ′ 2iπ lim K • γn in the sense of measures: [K(∂ t ) u] (t) = ⟨µ, ψ(t, .)⟩ , ( 10
) where ψ(t, ξ) is solution of the following evolution problem on (t, ξ) ∈ R * + ×R (of diffusive type):

∂ t ψ(t, ξ) = γ(ξ) ψ(t, ξ) + u(t), ψ(0, ξ) = 0. ( 11
) 1 Possibly at infinity 2 W 1,∞ loc (R; C) is the topological space of measurable functions f : R → C such that f, f ′ ∈ L ∞ loc (
→ ∫ t 0 u(s) ds, denoted ∂ -1 t . We indeed have (∂ -1 t u)(t) = ⟨δ, ψ(t, .)⟩ = ψ(t, 0), with ∂ t ψ(t, 0) = u, ψ(0, 0) = 0.
Beyond the measure framework, the general space of γsymbols is a quotient space of distributions, denoted ∆ ′ γ ; it is the topological dual of the space ∆ γ ∋ ψ(t, .) [START_REF] Montseny | Représentation diffusive[END_REF]). The composition product of operators has an equivalent in ∆ ′ γ , denoted ♯ γ or simply ♯: if µ and ν are respective γ-symbols of H(∂ t ) and 11) and µ the γ-symbol of K(∂ t ). In the same way, ∆ ′ γ can be extended to an algebra, denoted Σ γ , composed of the γ-symbols of operators γdiffusive of degree n ∈ N, the γ-symbols ν of K(∂ t ) • ∂ n t being characterized by the relation:

K(∂ t ), then µ♯ γ ν is a γ-symbol of H(∂ t ) • K(∂ t ).
µ = ν♯δ n (13)
where

δ n = δ♯δ♯...♯δ n times ∈ ∆ ′ γ is a γ-symbol of ∂ -n t .
The inversion of γ-symbols cannot be defined in ∆ ′ γ because this algebra is not unitary; this operation is nevertheless well-defined in Σ

γ . If µ ∈ Σ γ is a γ-symbol of K(∂ t ) such that K(∂ t ) -1 • ∂ -n t has a γ-symbol ν ∈ ∆ ′ γ , then ν = µ -1 ♯δ n and we have: [K(∂ t ) -1 u](t) = ⟨ µ -1 ♯δ n , ∂ n t ψ(t, .) ⟩ , ( 14
) with ψ solution of (11). Note in particular that operator ∂ -1 t defined above is the unique inverse of the derivative operator ∂ t . See [START_REF] Casenave | Représentation diffusive et inversion opératorielle pour l'analyse et la résolution de problèmes dynamiques non locaux[END_REF]) for details about the inversion of γ-symbols.

Summary

Let γ be a closed simple arc in C -verifying the property ( 7). An operator H(∂ t ) γ-diffusive of degree n admits the following state-realization:

∂ t ψ(t, ξ) = γ(ξ)ψ(t, ξ) + u(t), ψ(0, .) = 0, (15) (H(∂ t )u)(t) = < µ, ∂ n t ψ(t, .) > ∆ ′ γ ,∆γ , (16)
where

µ ∈ ∆ ′ γ is the γ-symbol of K(∂ t ) := H(∂ t ) • ∂ -n t ;
the essential conditions the operator has to satisfy to admit such a representation are:

• K holomorphic in Ω + γ , ( 17
) • K(p) → 0 when |p| → +∞ in Ω + γ , ( 18
)
K being the Laplace-symbol of K(∂ t ) given by:

K(p) = H(p) p n . ( 19
)
Thanks to the sector condition (7) verified by γ, the state realization is of diffusive type: so, cheap and precise numerical approximations of (15,16) can be easily built.

About numerical approximations

The state equation ( 11) is infinite-dimensional. To get numerical approximations, we consider M L a sequence of Ldimensional spaces of atomic measures on suitable meshes {ξ L l } l=1:L on the variable ξ; L-dimensional approximations µ L of the γ-symbol µ ∈ ∆ ′ γ are then defined in the sense of atomic measures, that is:

µ L = ∑ L l=1 µ L l δ ξ L l , µ L l ∈ C. ( 20
) If ∪ L M L is dense in the topological space ∆ ′ γ (that is, concretely, if ∪ L {ξ L l } is dense in R)
, then we can have [START_REF] Montseny | Représentation diffusive[END_REF]):

⟨ µ L , ψ ⟩ -→ L→+∞ ⟨µ, ψ⟩ ∀ψ ∈ ∆ γ ; (21) so, we have the following L-dimensional approximate state formulation of K(∂ t ) (with γ-symbol µ):    ∂ t ψ(t, ξ L l ) = γ(ξ L l ) ψ(t, ξ L l ) + u(t), l = 1 : L, ψ(0, ξ L l ) = 0 [K(∂ t ) u](t) ≃ ∑ L l=1 µ L l ψ(t, ξ L l ). (22)
One of the properties of the approach presented above is that most of non rational operators encountered in practice can be closely approximate with small L (see for example [START_REF] Montseny | Simple approach to approximation and dynamical realization of pseudodifferential time-operators such as fractional ones[END_REF])). In the context of identification of dynamic operators, this will be a great advantage because only a few numerical parameters µ L l will have to be identified from experimental data, while the property (21) will ensure the well-posedness and the robustness of the problem as soon as the operator to be identified admits a γ-symbol in Σ γ .

IDENTIFICATION OF AN OPERATOR BY MEANS OF ITS γ-SYMBOL

We consider here the problem of identification of an integral operator H(∂ t ) :

X -→ v = H(∂ t )X.
(23) Let v m and X m be some (direct or indirect and possibly noisy) measurements of the output v and the associated input X; the problem consists then in identifying the γsymbol of H(∂ t ) from these data.

Principle

Assume that H(∂ t ) is γ-diffusive of degree n and let µ be the γ-symbol of

K(∂ t ) = H(∂ t ) • ∂ -n t . We then have (see section 2): v = H(∂ t )X = ⟨µ, ∂ n t ψ X ⟩ = ⟨ µ, ψ ∂ n t X ⟩ .
(24) By denoting A ∂ n t X the operator defined by:

A ∂ n t X : µ -→ ⟨ µ, ψ ∂ n t X ⟩ , ( 25
)
we get a new formulation, linear with respect to µ: v = A ∂ n t X µ, (26) and which can also be written:

v m = A ∂ n t Xm µ + ϵ(µ), ( 27
)
where ϵ(µ) is the error equation. The least squares estimator μ of µ is then given by:

μ = arg min µ∈E ∥ϵ(µ)∥ 2 F , ( 28
)
where E is an Hilbert subspace of ∆ ′ γ and F another Hilbert space chosen a priori. The solution of ( 28) is formally obtained by:

μ = A † ∂ n t Xm v m , ( 29
)
where

A † ∂ n t Xm designates the pseudo-inverse of A ∂ n t Xm
(Ben-Israel and Greville [2003]).

In the sense of the hilbertian norm of F, the estimator μ is optimal. When the measurement X m is noisy, the estimator μ is also biased because A † ∂ n t Xm depends on the measurement noise. To mitigate this problem, it can be interesting to consider some classical bias reduction methods as the ones used in time-continuous system identification.

By denoting K(∂ t ) the operator of γ-symbol μ and

Ĥ(∂ t ) = K(∂ t ) • ∂ n t , the identified model is then written v = Ĥ(∂ t )X, or, under diffusive formulation: { ∂ t ψ = γ ψ + ∂ n t X, ψ(0, .) = 0 v = ⟨ μ, ψ ∂ n t X ⟩ . ( 30
)
Remark 2. Note that in the case where v m = v (no output measurement noise), we can get an unbiased estimation by identifying the operator H(∂ t ) -1 of input v and output X instead of H(∂ t ). We suppose

H(∂ t ) -1 is γ-diffusive of degree n ′ and we denote ν ∈ ∆ ′ γ the γ-symbol of K(∂ t ) = H(∂ t ) -1 • ∂ -n ′ t
. We then have:

X = H(∂ t ) -1 v = A ∂ n ′ t v ν, ( 31
)
and so get, after identification:

ν = A † ∂ n ′ t v X m . ( 32
)
The identified model then writes

v = K(∂ t ) -1 • ∂ -n ′ t X, that is, under diffusive formulation: { ∂ t ψ = γ ψ + ∂ n t X, ψ(0, .) = 0 v = ⟨ (ν) -1 ♯δ n ′ +n , ψ ∂ n t X ⟩ . ( 33
)
The computation of (ν) -1 ♯δ n ′ +n can be numerically performed as shown in [START_REF] Casenave | Représentation diffusive et inversion opératorielle pour l'analyse et la résolution de problèmes dynamiques non locaux[END_REF]).

More information about this identification method will be found in [START_REF] Casenave | Représentation diffusive et inversion opératorielle pour l'analyse et la résolution de problèmes dynamiques non locaux[END_REF]Montseny [2009b, 2010]).

Prefiltering with an invertible convolution operator

The identification model ( 26) can be equivalently transformed by composition with any invertible causal convolution operator Q(∂ t ):

Q(∂ t )v = Q(∂ t ) • H(∂ t )X = H(∂ t ) • Q(∂ t )X. ( 34
)
By denoting ṽ = Q(∂ t )v and X = Q(∂ t )X, the model is then written: ṽ

= H(∂ t ) X. ( 35
)
When applying the identification method to model ( 35), the estimator of µ is written:

μ = A † ∂ n t Xm ṽm , ( 36
) with Xm = Q(∂ t )X m and ṽm = Q(∂ t )v m .
When n ̸ = 0, this transformation is necessary to make the identification process possible, otherwise, the amplification of the noise measurement of the term ∂ n t X m in ( 28) and ( 29) would make the identification impossible. So, the operator Q(∂ t ) has to be chosen so that the measurement noise is not amplified by the term ∂

Q(p) = σ n (p + σ) n , ( 37
)
where σ > 0 (the cutoff frequency) will be chosen in such a way that ∥ϵ(µ)∥ 2 F is as "small" as possible. Remark 3. The transfer function Q(p) could also be optimized in order to minimize the estimation error.

Case of multiple trajectories

Consider a set of input trajectories X j , j = 1 : J and the associated output trajectories v j = H(∂ t )X j . Let v j m and X j m be some measurements of v j and X j . Then, without any change of notations, model ( 26) can be extended to the general case of multiple trajectories simply by defining:

v = (v 1 , ..., v J ) T , X = (X 1 , ..., X J ) T , v m = (v 1 m , ..., v J m ) T , X m = (X 1 m , ..., X J m ) T ,
and

A ∂ n t X : µ -→    < µ, ψ ∂ n t X 1 > . . . < µ, ψ ∂ n t X J >    .
(38)

Recursive formulations

Recursive formulations of (29) can be established under the form [START_REF] Garcia | Identification of the dynamics of a lead acid battery by a diffusive model[END_REF]):

μt = μt-∆t + K t-∆t (v m -A ∂ n t Xm μt-∆t ) |[0,t] ; (39)
such formulations allow real-time identification (or even the pursuit of µ in case of slowly varying operators H(t, ∂ t )).

In the numerical point of view

Let (v, X) be such that v = H(∂ t )X. We consider a discrete data set {v k m , X k m } k=0:K where v k m and X k m are the respective measurements of v and X at time instant t k . We also denote v m and X m the continuous measured trajectories such that v m (t k ) = v k m and X m (t k ) = X k m . The time instants t k are given by:

t 0 = 0, t k = t k-1 + ∆t k , k = 1 : K. ( 40
)
We consider a discretization {ξ L l } l=1:L of the variable ξ and the associated L-dimensional approximation µ L of µ given by (20). At time instant t k , we consider the following identification model of finite dimension:

v k m = φ T m (t k ) µ + ϵ(t k , µ), ( 41
) where :

• µ = (µ L 1 , ..., µ L L ) T , • φ T m (t k ) = [Ψ k 1 , ..., Ψ k L ] with Ψ k l = ψ ∂ n t Xm (t k , ξ L l ), • and ϵ(t k , µ) is the error equation at time t k .
The least squares estimator μ is then given by: μ = arg min

µ∈C L K ∑ k=1 ϵ 2 (t k , µ)∆t k , ( 42
)
that is:

μ = [ K ∑ k=1 φ m (t k ) φ T m (t k )∆t k ] -1 K ∑ k=1 φ m (t k ) v k m ∆t k . (43)

On the choice of γ

To choose the contour γ on which depends the problem (28), recall first that γ has to verify the sector condition (7), and then that the transfer function

K of K(∂ t ) = H(∂ t ) • ∂ -n t
has to be analytic in Ω + γ (see( 17)). So all the singularities of K have to be inside the domain Ω - γ delimited by γ. However, as the operator H(∂ t ) is unknown, so is the position of the singularities of K. As a consequence the contour γ will be chosen in such a way that the domain Ω - γ is sufficiently big to contain all the singularities of K. In practice, we often take a contour of sector type: γ(ξ) = |ξ| e i sign(ξ)( π 2 +α) , (44) with α ∈]0, π 2 ] a sufficiently "small" angle. Note that, however, the more γ is close to the axis iR (stability limit axis), the finer (and so the more expensive numerically) the discretization in ξ has to be in order to get a good approximation of the model. So in practice, we first identify the system with a small value of angle α and a great value L. Then, if the identification results are good with these values, we can iterate the process with a greater angle α (α π 2 ) and a smaller L. Note that in practice, the available information about the operator H(∂ t ) can help us in the choice of γ.

Some indications about the choice of the discretization points ξ L

l can be found in [START_REF] Montseny | Représentation diffusive[END_REF]. Note however that, from the discrete data {v k m , X k m }, we can identify the impulse response H(iω) only in the frequency band

[ 2π t K ; 2π 2∆t ],
where ∆t = max(∆t k ). Consequently, the band [ξ L 1 ; ξ L L ] covered by the discretization in ξ will be chosen in such a way that5

[ 2π tK ; 2π 2∆t ] ⊂ γ([ξ L 1 ; ξ L L ]) .

IDENTIFICATION OF VOLTERRA MODELS

Now we consider a nonlinear Volterra model of the form:

H(∂ t )X = f (u, X) + v, ( 45
) where H(∂ t ) is a causal convolution operator and f is a nonlinear function. The problem under consideration consists in identifying both operator H(∂ t ) and function f from some measurements u m , v m and X m of u, v and X. We will describe two identification methods: with the first one, H(∂ t ) and f are identified simultaneously, whereas with the second one, H(∂ t ) is first identified alone after cancellation of the non linear term f in the model. Only the theoretical principles of these two methods are described; the numerical implementation is similar to the one described in section 3.5.

Suppose H(∂ t ) is γ-diffusive of degree n and let µ be the

γ symbol of K(∂ t ) = H(∂ t ) • ∂ -n t .
Then, (45) can be rewritten:

A ∂ n t X µ = f (u, X) + v. ( 46
)

Simultaneous identification of H(∂ t ) and f

We consider a topological basis {g p ⊗ k q } p,q=1:+∞ of a tensorial product of Hilbert spaces to which f belongs. We so have:

f = ∑ p,q g p ⊗ k q a pq , with a := (a pq ) p,q ∈ ℓ 2 . ( 47
)
The Volterra model ( 45) can then equivalently be expressed under the linear form:

G u,X (µ, a) = v, ( 48
) with G u,X the operator defined by:

G u,X : (µ, a) → A ∂ n t X µ - ∑ p,q g p (u) k q (X) a pq . ( 49
)
The least squares estimator (μ, â) of (µ, a) is defined by:

(μ, â) = arg min (µ,a)∈E×ℓ 2 ∥ϵ(µ, a)∥ 2 F , ( 50
)
where ϵ(µ, a) = G um,Xm (µ, a) -v m is the error equation. Formally we so get:

(μ, â) = G † um,Xm v m . ( 51
)
When n ̸ = 0, a prefiltering with a convolution operator (see section 3.2) is necessary so that the noise is not amplified.

For complementary details about this method, we will refer to [START_REF] Casenave | Identification of Nonlinear Volterra Models by means of Diffusive Representation[END_REF]).

Identification of H(∂ t ) after cancellation of f

For simplicity, we here suppose that f (u, X) = f (X).

Let also assume that n < 2 (note that in other cases, the identification method presented here-after can not be applied due to the amplification of the measurement noise6 by operator ∂ n t : a prefiltering with a convolution operator cannot mitigate the problem because the method is based on the fact that f is a function, which is not the case of

Q(∂ t ) • f ). Let ε ∈ R + and (x, y) ∈ ( C 0 ([0, T ]; X) ) n × ( C 0 ([0, T ]; Y) ) n
with X, Y two Banach spaces; we denote:

Ω x,ε := ∪ i,j=1:n {(i, j)} × Ω i,j x,ε , (52)
where Ω i,j x,ε :=

{ (t, τ ) ∈ [0, T ] 2 ; x i (t) -x j (τ ) ε } (see figure 1 for some examples). Fig. 1. Examples of (t, τ ) in Ω i,i
x,ε (top) and in Ω i,j

x,ε

(bottom) Definition 5. The ε-cancellation operator D ε is defined by: (x, y) → D ε (x, y) : Ω x,ε → Y (i, j, t, τ ) → y i (t) -y j (τ ). ( 53
)
Proposition 4. (1) D ε is continuous.

(2) For any continuous function f : X 0 ⊂ X → Y, we have the cancellation property:

D 0 (x, f • x) = 0. ( 54
)
Thanks to proposition 4 and by application of operator D 0 (X, .) to both members of the equation, model ( 45) is rewritten:

{ ⟨ µ, D 0 ( X, ψ ∂ n t X )⟩ = D 0 (X, v) ⟨ µ, ψ ∂ n t X (0, .) ⟩ -f (X(0)) = v(0). ( 55
)
The interest of this new formulation is that, on the one hand, up to the quantity f 0 := f (X(0)), the model is independent of f , and on the other hand, both µ and f 0 are linearly involved. Formally, these unknowns can then be determined by means of the pseudo-inverse of operator:

Y 0 : (µ, f 0 ) → ( ⟨ µ, D 0 (X m , ψ ∂ n t Xm ) ⟩ ⟨ µ, ψ ∂ n t Xm (0, .) ⟩ -f 0 ) . ( 56
)
In general however, on the one hand the set Ω Xm,0 may be too poor to get the strict equivalence of ( 45) and ( 55) and, on the other hand, the available data X m and v m are noisy. In practice, we so consider the weakened problem:

min (µ,f0) ∥ϵ(µ, f 0 )∥ 2 , ( 57
)
where the equation error ϵ(µ, f 0 ) is given by:

ϵ(µ, f 0 ) = Y ϵ (µ, f 0 ) - ( D ε (X m , v m ) v(0) ) , (58)
with:

Y ϵ : (µ, f 0 ) → ( ⟨ µ, D ε (X m , ψ ∂ n t Xm ) ⟩ ⟨ µ, ψ ∂ n t Xm (0, .) ⟩ -f 0 ) . ( 59
)
The estimator (μ, f0 ) of (µ, f 0 ) is then defined by:

(μ, f0 ) = Y † ε ( D ε (X m , v m ) v m (0) ) . ( 60
)
By means of standard regression methods, the function f can then be easily estimated from its "pseudo graph":

P f = ∪ k {( X k m , ⟨ μ, ψ ∂ n t Xm (t k , .) ⟩ -v k m )} . ( 61
)
The hypothesis made on f and n are more restrictive than those made for the method of simultaneous identification of f and H(∂t). Nevertheless, this method presents some interesting advantages. First, as the identification of f is uncoupled from the one of H(∂t), it is possible to adapt the choice of basis functions to the behavior of function f given by its pseudo-graph. Second, this method can also be used when the function f changes from one data set to the other. This identification method is still under study: more information about it will be found in (Casenave and Montseny [2009a]).

Iterative method

The identified quantities μ and f can be used as initial values for an iterative identification method in order to improve the estimation quality. For example, let's assume that the estimated operator Ĥ(∂ t ) is exact; we can identify the function f again, this time independently from H(∂ t ).

From the so-identified function denoted f1 , we identify H(∂ t ), and we then reiterate until the quality of the estimation becomes sufficient. Note that the convergence of such a method is still under study.

Identified model

The physical system under consideration can be described by the model Ĥ(∂ t )X = f (u, X) + v where Ĥ and f are deduced from the identified parameters. Furthermore, by assuming that Ĥ(∂ t ) -1 is γ-diffusive of degree n ′ , we deduce from the equivalent expression

X = Ĥ(∂ t ) -1 ( f (u, X) + v ) the following input-output state realization :    ∂ t ψ = γ ψ + f ( u, ⟨ (μ) -1 ♯δ n+n ′ , ψ ⟩) + v, ψ(0, .) = 0, X = ⟨ (μ) -1 ♯δ n+n ′ , ∂ n ′ t ψ ⟩ , ( 62
) which is available up to the computation of (μ) -1 ♯δ n+n ′ which can be numerically performed as shown in [START_REF] Casenave | Représentation diffusive et inversion opératorielle pour l'analyse et la résolution de problèmes dynamiques non locaux[END_REF]).

APPLICATION TO A MODEL OF FLAME

In this section, we illustrate the identification methods presented above by implementing them on data elaborated from numerical simulations of a complex dynamic phenomenon [START_REF] Joulin | Point-source initiation of lean spherical flames of light reactants: an asymptotic theory[END_REF]).

The model under consideration

In [START_REF] Joulin | Point-source initiation of lean spherical flames of light reactants: an asymptotic theory[END_REF]), Joulin elaborated a Volterra model to describe, in suitable thermodynamic conditions, the evolution of a spherical flame initiated by a source at point 0 in a mixture of reactive species. Under some reasonable physical hypothesis, such a phenomenon can be described by a system of two partial differential equations relating to the temperature and the mass density of the mixture. By considering the reactive zone as a thin sheet located on a sphere with radius X(t), Joulin has established that when the flame is developing in free space, X is solution of the following nonlinear singular Abel-Volterra equation 7 (e(t) designates the source strength at time t):

X(t)∂ 1 2 t = 2 X(t) ln X(t) + 2 e(t) ∀t > 0, ( 63
) with the additional conditions: X(0) = 0, e 0, X 0 (whose physical interpretation is obvious). By denoting

H(∂ t ) = 1 2 ∂ 1 2
t , f (u, X) := ln X = f (X) and u = e X , (63) can be formally rewritten under the form (45). We have:

• H(p) p = 1 2 p -1 2 holomorphic in C R -,
• and H(p) p → 0 when |p| → 0 in C R -; so H(∂ t ) is γ-diffusive of degree 1 for any γ of the sector form (44).

As shown in [START_REF] Audounet | A threshold phenomenon in the propagation of a point source initiated flame[END_REF]), the evolution problem ( 63) is well-posed, that is the solution X exists, is unique and depends continuously on e. It has also been shown that there exists a threshold relating to the power of the source e, beyond which the flame is developing whereas a quenching occurs below, as it is highlighted in Fig 2a and2b, the source function used for the simulation being given as in [START_REF] Audounet | A threshold phenomenon in the propagation of a point source initiated flame[END_REF] by:

e(t) = E t 0.3 (1 -t) 1 [0,1] (t). (64)
In real conditions, various perturbations are involved in the evolution of X (due for example to the loss of spatial symmetries), and both the convolution operator H(∂ t ) and the function f will be more or less far from the ideal ones. So, an identification process can be justified if accuracy of the model is required.

Numerical identification results and comments

We consider in the sequel the problem of identification of H(∂ t ) and/or f from data (e m , X m ) obtained from highly accurate numerical simulations of (63) [START_REF] Montseny | Desingularization of non local dynamic models by means of operatorial transformations and application to a spherical flame model[END_REF]).

The data are composed of 9 discrete simulated trajectories (e j,k m , X j,k m ) k=0:K , j = 1 : 9 (see figure 2c), associated with 9 different sources e j of the form8 (64) with E = 1.5, 1.6, 1.7, 1.73, 1.74, 1.75, 2.0, 3.0 and 5.0 respectively. The time step ∆t is constant and equals 10 -5 and the final time is given by t K = 6 (K = 6 × 10 5 ). The inputs are supposed to be known: e j,k m = e j (t k ), (65) and the outputs are corrupted by additive measurement noise:

X j,k m = X j (t k ) + v j,k , (66) where {v j,k } is a numerical white noise with zero mean and of standard deviation s = 10 -1 .

To identify H(∂ t ), we use L = 80 values ξ L l geometrically spaced to cover 5 decades from 10 0 to 10 5 , and we consider a contour γ of the form (44) with α = π 4 . To identify f , we use 30 basis functions k q , q = 1 : 30, of hat type and one function basis g 1 = 1.

The frequency response of the identified operator Ĥ(∂ t ) and the curve of the identified function f are given in figure 3 when f is supposed to be known (see section 3), when f and H(∂ t ) are simultaneously identified (see section 4.1) and when H(∂ t ) is first identified alone after cancellation of f (see section 4.2). For the 3 methods, the identification of H(∂ t ) is good in the frequency band [10 0 , 10 2 ] ⊂ [ 2 π K ∆t , 2 π 2 ∆t ] = [1.05, 3.14 10 5 ] (in the free noise case, this frequency band covers 4 decades). The identification of f is also good, even around the singularity of the function f at X = 0. Note however that the greater the value of X is, the larger the difference between f (X) and its estimate f (X) is; this is partly due to the fact that there are more data X j,k m with value near from zero than others.

CONCLUSION

In this paper, we have shown how to use the diffusive representation in order to identify the dynamic part of a system. Thanks of the linearity of the γ-symbol in the state realization of the operator H(∂ t ), we can identify it by means of classical least squares methods. Several questions must yet be studied in order to improve the results. For example, among the most significant, the involved hilbertian norms should be judiciously chosen and adapted to the specific properties of the class of models under consideration; indeed, this choice is crucial in terms of sensitivity with respect to perturbations of any nature. The estimation bias could also be reduced. All these questions are currently under study.

  Source (---) and radius of the flame (-) with E = 1.7390.

  Source (---) and radius of the flame (-) with E = 1.7393.

  Trajectories X used for identification.

Fig

  Fig. 2.

Fig. 3 .

 3 Fig. 3. Identification results.

  that is f and f ′ are locally bounded in the almost everywhere sense).3 This convergence mode means that on any bounded set P , γn |P -

	Definition 4. The measure µ defined in theorem 3 is called
	γ-symbol of operator K(∂ t ). The function ψ solution of
	(11) is called the γ-representation of u.
	Remark 1. Note in particular that thanks to (10), the
	Dirac measure δ is clearly a γ-symbol of the operator
	u
	γ |P → 0 and γ′ n |P -γ ′ |P → 0 uniformly.

  Basically, it behaves like a nth order low-pass filter; here we simply consider the filter Q(∂ t ) with transfer function:

	n t Xm . It has to suf-ficiently attenuate high frequencies (i.e. |Q(iω)| ∼ H.F 1 ω n ),
	without amplifying low and middle ones (i.e. |Q(iω)| ∼ L.F
	1).

In the sense of the strong topology of ∆ ′ γ .

The quantity |γ(ξ)| is comparable to a frequency. In the case of a contour of sector type we have |γ(ξ)| = |ξ|.

In the case where n = 1, we have ψ ∂tX = ∂tψ X = γψ X + X; so the measurement noise is not amplified by the computation of ψ ∂tX .

Note that the source function is physically realistic but rather poor from the point of view of information, which strengthen the difficulty of the problem.


Céline Casenave
email: casenave@laas.fr

Gérard Montseny
email: montseny@laas.fr

Henri Camon
email: camon@laas.fr

François Blard
email: francois.blard@laas.fr

François Blard Identification

Identification of dynamic nonlinear thermal transfers for precise correction of bias induced by temperature variations

à 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 Dynamic Nonlinear Thermal Transfers for Precise Correction of Bias induced by

Temperature Variations Céline Casenave * † , Gérard Montseny * † , Henri Camon * † and Franc ¸ois Blard * †

I. INTRODUCTION

In general, high precision devices, such as MEMS devoted to metrology applications for example, need to include algorithms allowing to correct biases resulting from variations of the ambient temperature. Most of time, static corrections are involved, which means that the ambient temperature (measured by means of a sensor) is supposed to be also the temperature at any point of the device itself. This is of course not exactly true because the ambient temperature is never strictly constant and so, due to the complex distributed nature of physical objects, the temperature inside the device evolves following a complex diffusion equation. So, when measuring any physical quantity (for example a voltage) on the device under consideration, small residual fluctuations due to temperature variations remain in spite of the static correction, because dynamic thermal phenomena are involved, which cannot be described by means of a static correspondence.

In some cases, when very high precision is required, it can be judicious to envisage to build some dynamic corrections, elaborated from a suitable dynamic treatment of the temperature data, in such a way that the macroscopic effects resulting from the distributed evolution of the temperature inside the device can be taken into account, and then corrected.

In this paper, we present a method devoted to such a dynamic correction.

II. RECALLS ON STANDARD (STATIC) POLYNOMIAL

CORRECTION

Consider a quantity U (to be measured), linked to the ambient temperature T through a static relation of the form:

U = f (T ). ( 1
)
The aim is to correct the effect of the temperature T , that is to get an (if possible optimal) estimation of the quantity

U 0 = f (T 0 )
, where T 0 is chosen a priori as the reference temperature. If f is regular (which is the case in concrete situations), U can be expressed following the Taylor development around T 0 :

U = f (T 0 ) + f (T 0 )(T -T 0 ) + 1 2 f (T 0 )(T -T 0 ) 2 + ..., (2)
that is, if we consider that terms beyond the n th order are negligible:

U = U 0 + n i=1 a i (T -T 0 ) i . ( 3
)
So if the coefficients a i have been previously accurately estimated (from experimental measurements) and T is known, the bias induced by sufficiently slow temperature variations around T 0 can be corrected thanks to the relation:

U 0 = U - n i=1 a i (T -T 0 ) i . ( 4
)
Indeed, when U is measured N times with some additive zero mean noise w, we then have Ũk = U k + w k , k = 1 : N and so, we get the correction relation:

Ũ0,k = Ũk - n i=1 a i (T k -T 0 ) i = U k - n i=1 a i (T k -T 0 ) i + w k (5) = U 0 + w k , k = 1 : N.
An estimate of the unknown quantity U 0 can then be obtained by filtering the data Ũ0k , for example by:

Û0 = 1 N N k=1 Ũk - n i=1 a i (T k -T 0 ) i . ( 6
)
We then have from [START_REF] Casenave | Identification of a Nonlinear Dynamic Model of MEMS with Unstable Switching Zone[END_REF] and thanks to the law of large numbers:

E( Û0 ) = U 0 , (7)
Û0 = 1 N N k=1 (U 0 + w k ) → U 0 when N → ∞; (8)
so, the estimate Û0 is convergent and unbiased: the effects of temperature variations have been corrected. Finally, if the unknown quantity U 0 is not exactly constant but presents some very slow (and small) variations, it can be pursued for example by means of the following moving average estimate based on the corrected data Ũ0,k :

Ū0,k = 1 N N -1 q=0 Ũk-q - n i=1 a i (T k-q -T 0 ) i , (9)
with the following residual estimation noise (inherited from the measurement noise):

ε k = 1 N N -1 q=0 w k-q . ( 10
)

Preliminary estimation of the correction coefficients from experimental data

The correction coefficients a i are estimated from experimental data, whose number must be sufficiently large to make the residual estimation error negligible.

In practice, this step can be very long, particularly if measurement data are significantly corrupted by noise. The classical estimation of coefficients a i from data Ũk , k = 1 : K is in general obtained from the least squares method, that is by solving the problem:

min ai∈R K k=1 Ũk -a 0 - n i=1 a i (T k -T 0 ) i 2 ; ( 11
)
this problem is rewritten under the matrix form:

min a∈R n+1 (Qa -Ũ ) T (Qa -Ũ ) (12)
and the optimal solution a * = (a * i ) is classically given by:

a * = (Q T Q + εI) -1 Q T Ũ (13)
with ε 0 a small parameter devoted to numerical conditioning of the matrix inversion.

III. PRECISE CORRECTION OF BIAS FROM IDENTIFIED DYNAMIC THERMAL TRANSFER

A. Dynamic influence of the ambient temperature on a physical quantity

In practice, a physical system such as a MEMS with its environment is a distributed system: when the ambient temperature T is not constant, the temperature of the system cannot be a simple scalar: it becomes a scalar field θ(t, x), where each (vector) value of the spatial variable x is associated with a particular physical point of the system. If the variations of T remain small, the evolution of the field θ is theoretically governed by a complex but linear partial differential equation with input T (t), of the abstract form:

∂θ(t, x) ∂t = A(∇) θ(t, x) + B(x) T (t), x ∈ Ω ⊂ R 3 , ( 14
)
where A(∇) is a second order linear differential operator, associated with suitable boundary conditions, and B is the input operator describing how the exterior temperature T acts on the evolution of θ. Because equation ( 14) is in general defined on t > t 0 , the initial condition must be added in order to get a well-posed problem:

θ(t 0 , x) = θ 0 (x), x ∈ Ω. (15)
Except in some ideal simple cases, such an equation is in practice impossible to formulate explicitly, due to the complexity of any physical system.

Similarly, a quantity U relating to the physical system under consideration and depending at first sight on T only, depends in general on the whole field θ. As in the ideal static cases, this dependence can of course be nonlinear, in such a way that we have formally:

U (t) = F(θ(t, x)), ( 16
)
where F is a nonlinear spatial operator. So, the quantity U appears as an output of ( 14) and the dynamic (nonlinear) transfer T -→ U is then resulting from the input-output correspondence defined by (14,16).

In most of cases, the evolution of T (t) is slow and it can be considered that the difference |θ(t, x)-T (t)| can be neglected and then, the static approach previously described can be used. In some cases however, namely when very high precision is required, small differences between θ(t, x) and T (t), which result from the dynamic nature of equation ( 14), can generate some significant residual biases which cannot be suppressed by static correction. In such cases, the static model (3) is no more sufficiently accurate and a dynamic correction is needed.

B. Dynamic correction of bias resulting from temperature variations

This correction method is based on an explicit and universal differential model of the dynamic transfer (14,16), as described here-after1 . First, thanks to the linear and diffusive nature of equation ( 14), a generic 1D diffusive input-state equation is used to model the underlying thermal dynamics. By assuming for simplicity that the mean ambient temperature is T 0 = 0, this model is:

⎧ ⎨ ⎩ ∂ψ(t, ξ) ∂t = -ξ ψ(t, ξ) + T (t), ξ > 0, t > t 0 ψ(t 0 , ξ) = 0. (17)
On the other hand, from linearity, the contribution of the initial condition θ 0 is separately expressed from the exponential family:

ψ 0 (t, ξ) = e -ξt ψ 0 (ξ); (18)
note that the function (t, ξ) → ψ(t, ξ) + ψ 0 (t, ξ) is solution of (17) with initial condition ψ 0 (ξ) in place of 0. We then define the function:

Ψ(t, ξ) := ψ(t, ξ) + ψ 0 (t, ξ) -T (t) (19)
which expresses the difference between the ambient "static" temperature and the dynamic field ψ(t, ξ) + ψ 0 (t, ξ); this function will be devoted to dynamic correction. From ψ (the solution of ( 17)) and ψ 0 , we can then model the unknown dynamic transfer T → U under the form U = G(T, Ψ) with G a suitable nonlinear operator; a generic formulation of G is the following:

U (t) = U 0 + i a i T (t) i + μ 1 (ξ) Ψ(t, ξ) dξ + μ 2 (ξ 1 , ξ 2 ) Ψ(t, ξ 1 ) Ψ(t, ξ 2 ) dξ 1 dξ 2 (20) + μ 3 (ξ 1 , ξ 2 , ξ 3 ) Ψ(t, ξ 1 ) Ψ(t, ξ 2 ) Ψ(t, ξ 3 ) dξ 1 dξ 2 dξ 3 + ...,
where the coefficient a i and the functions μ i allow to synthesize from T (t) and the field Ψ(t, ξ) a wide class of nonlinear transfers. Thanks to the nice properties of equation ( 17), an approximated synthesis with only a few Ψ j (t) := Ψ(t, ξ j ), ξ j > 0, of the form:

U (t) = U 0 + n i=1 a i T (t) i + J1 j=1 μ 1,j Ψ j (t) + J2 j=1 J 2 j =1 μ 2,j,j Ψ j (t) Ψ j (t) + J3 j=1 J 3 j =1 J 3 j =1 μ 3,j,j Ψ j (t) Ψ j (t) Ψ j (t) + j j j j ..., (21)
will be in practice sufficient to get accurate approximations of the dynamic nonlinear transfer T → U . Finally, if the temperature variations are sufficiently small (this will be the case in the problem studied in the next section), the cross terms Ψ j Ψ j , Ψ j Ψ j Ψ j , etc., can be neglected and the following simplified synthesis will be yet sufficient to get significant improvement compared to the static model (3):

U (t) = U 0 + n i=1 a i T (t) i + J j=1 μ j (ψ j (t) + ψ 0,j (t) -T (t)), (22)
in such a way that (22) can be rewritten, from (18) and with ν j := μ j ψ 0 (ξ j ), under the following form in which all the correction parameters a i , μ j , ν j are linearly involved:

U (t) = U 0 + n i=1 a i T (t) i + J j=1 μ j (ψ j (t)-T (t))+ J j=1
ν j e -ξj t .

(23) Algorithm for dynamic correction Similarly to the static case, when U is measured at times t k = kΔt with some additive zero mean noise w k that is Ũk = U k +w k , and thanks to the fact that e -ξj t k → 0 when k → +∞, we finally get the asymptotic dynamic correction formula, suitable for large k (NB: the correction coefficients are a i , μ j ):

Ũ0,k = Ũk - n i=1 a i (T k ) i - J j=1 μ j (ψ j,k -T k ), (24)
where the quantities ψ j,k , on which is based the dynamic correction, are computed from the temperature data T k via the following discrete-time dynamic relation deduced from integration of (17):

ψ j,k+1 = e -ξj Δt ψ j,k + 1-e -ξ j Δt ξj T k . ( 25
)
Finally, exactly as in the static case, Ũ0,k = U 0 + w k and so, estimates Û0 or Ū0,k of the quantity U 0 can yet be taken as:

Û0 = 1 N N k=1 Ũ0,k , Ū0,k = 1 N N -1 q=0 Ũ0,k-q ; (26)
we have yet, from (24) and thanks to the law of large numbers:

E( Û0 ) = U 0 , ( 27
) Û0 = 1 N N k=1 (U 0 + w k ) → U 0 when N → ∞, ( 28
)
ε k = 1 N N -1 q=0 w k-q . ( 29
)

Preliminary estimation of the correction coefficients from experimental data

The correction coefficients a i , μ j , ν j are estimated from experimental data whose number must be sufficiently great to get negligible estimation errors. This estimation from data (with noise) Ũk , k = 1 : K is obtained again by means of the least squares method, by minimizing with respect to a i , μ j , ν j ∈ R the quantity:

J = K k=1 Ũk -a 0 - n i=1 a i (T k ) i - J j=1 μ j (ψ j,k -T k ) - J j=1 ν j e -ξt k ⎞ ⎠ 2 , ( 30
)
which can be rewritten under matrix form, with

(a, μ, ν) ∈ R (n+1) × R J × R J : min (a,μ,ν) ⎡ ⎣ Q ⎛ ⎝ a μ ν ⎞ ⎠ - Ũ ⎤ ⎦ T • ⎡ ⎣ Q ⎛ ⎝ a μ ν ⎞ ⎠ - Ũ ⎤ ⎦ , ( 31
)
the solution (a * , μ * , ν * ) of which is given by:

⎛ ⎝ a * μ * ν * ⎞ ⎠ = Q T Q + εI -1 Q T Ũ. (32)

IV. APPLICATION TO A HIGH PRECISION VOLTAGE REFERENCE BASED ON AN ELECTROSTATICALLY ACTUATED MEMS

In the sequel, we present some experimental results obtained by use of the method described above. The physical system under consideration is an electrostatically actuated MEMS devoted to the construction of a voltage reference of high precision and stability. From the electrical point of view, it is simply a variable capacitance in which the electrostatic force is opposed to the one of a mechanical spring between the two electrodes. So, there exists a so-called pull-in voltage, beyond which the mechanical force is no more able to balance the electrostatic one: this defines the reference voltage, which is determined by the only mechanical design of the MEMS [START_REF] Suhonen | AC and DC voltage standards based on silicon micromechanics[END_REF], [START_REF] Kärkkäinen | MEMS based AC voltage reference[END_REF], [START_REF] Kyynäräinen | Stability of microelectromechanical devices for electrical metrology[END_REF].

For stability tests, the electrical environment of the MEMS is shown (in a simplified form) in Fig. 2. The DC voltage applied between the MEMS electrodes is chosen slightly less than the pull-in voltage. Variations of the ambient temperature generate variations of the electrical capacitance of the MEMS (around 10pF) and therefore variations of the measured AC voltage.

In order to evaluate the long time stability of the electrical characteristics of the device and its electrical environment, it is necessary to correct these voltage variations. However, due to the particular geometry of the MEMS whose desired electromechanical properties result from a complex design with very thin spring suspension of the mobile electrode (see Fig. 1) and also to the imperfections of the other associated devices (see Fig. 2), we can expect that large time constants are possibly involved in the influence of temperature on the AC voltage measurement. In such a case, the temperature cannot be uniform in the whole system if the ambient temperature variations are not sufficiently slow. Taking into account that very high precision is required, the dynamic correction presented above has therefore been implemented for a better correction. The efficiency of the method is illustrated by comparing the results with those obtained from standard static correction.

The measured data are shown in Fig. 3; the sample period is Δt = 15s. Significant variations of the AC voltages are visible; they mainly result from ambient temperature variations. For static correction, we have taken n = 1 (linear correction) or n = 2. For dynamic correction, J = 6 time constants 1/ξ j have been used.

In Fig. 4, we can see both the measured data and a dynamic correction a 0 -

n i=1 a i (T k ) i - J j=1 μ j (ψ j,k -T k )
, where a 0 is the estimate of the exact (unknown) AC voltage U 0 . The correction is then subtracted to the measured data to get the socalled residual variations. These residual variations are shown in Fig. 5,6, either non smoothed, or smoothed by means of a standard moving average filter. We can clearly see that a significant improvement is obtained from dynamic correction.

Evolutions of the temperature T and the "auxiliary temperatures" ψ j +ψ 0,j are shown in Fig. 7,8,9. In Fig. 10, we can see the contribution to these evolutions of the terms inherited from initial conditions, that is: ψ 0,j (t k ) = νj μj e -ξj t k . Finally, Fig. 11 shows the differences between the auxiliary temperatures and T , that is the functions Ψ j = ψ j + ψ 0,j -T .

From a more quantitative point of view, we have the following results:

V. CONCLUSION

As clearly highlighted by these quantitative results, dynamic correction can in some cases generate significant improvements. Some residual variations remain after correction; they probably cannot be entirely induced by noises only. Indeed, great deviations of Ψ j (t) sometimes appear (see Fig. 11), which suggest that the only linear part of the dynamic correction is possibly no more sufficient in such cases. So, it could be judicious to involve higher order (nonlinear) terms of the series (21), in particular quadratic ones: μ 2,j,j Ψ j (t) Ψ j (t). This will be studied in a further work, with a new set of data devoted specifically to this problem.

  Mean measured AC voltage: 0.1040mV, Standard deviation before correction of the AC voltage: 1.114 × 10 -3 mV.Standard deviation of the AC voltage after correction (mV):non smoothed smoothed stat. lin.9.985 × 10 -5 4.911 × 10 -5 stat. nonlin. 9.780 × 10 -5 4.471 × 10 -5 dyn. lin.6.769 × 10 -5 1.331 × 10 -5 dyn. nonlin. 6.725 × 10 -5 1.173 × 10 -5 .Reduction of the residual variations of the AC voltage (smoothed corrected data):stat. lin.-27, 3 dB stat. nonlin. -27.9 dB dyn. lin.-38.6 dB dyn. nonlin. -39.5 dB.Improvement by use of dynamic correction versus static one: linear -11.3 dB nonlinear -11, 6 dB.

Fig. 1 .Fig. 3 .

 13 Fig. 1. The electrostatically actuated MEMS under consideration

Fig. 4 .Fig. 5 .Fig. 6 .Fig. 7 .Fig. 8 .Fig. 9 .Fig. 10 .Fig. 11 .

 4567891011 Fig. 4. Measured AC voltage (:), dynamic correction (-)

This approach approach is relating to the so-called diffusive representation. More details about this theory and its applications in identification of nonlinear complex models or other various fields can be found in [3],[START_REF] Montseny | Représentation diffusive[END_REF].


C Casenave

G Montseny

OPTIMAL IDENTIFICATION OF DELAY-DIFFUSIVE OPERATORS AND APPLICATION TO THE ACOUSTIC IMPEDANCE OF ABSORBENT MATERIALS

Keywords: Optimal identification, diffusive operators, diffusive representation, acoustic impedance, absorbent wall

come

Optimal identification of delay-diffusive operators and application to the acoustic impedance of absorbent materials

Céline Casenave,

INTRODUCTION

Identification of dynamic input-output systems is a central problem each time the knowledge of the process under consideration is too imprecise to get sufficiently accurate models from physical analysis (note anyway that such an analysis is sometimes impossible to perform because too complex). However, identification problems can also present some serious difficulties. Namely, when the input-output operator under consideration involves distributed underlying phenomena, significant non rational components are quasi systematically generated, which are associated, in the time domain, to infinitedimensional realizations.

When the operator to be identified is linear, a convenient and rather general approach consists in working in the frequency domain. So, any causal operator can be well-defined by its symbol Q (t, iω), that is the Fourier transform of the so-called impulse response h(t, s) of the operator. In the convolutive case, the symbol no more depends on t and the problem of identifying Q(iω) can be classically solved from physical measurements by means of Fourier techniques. Note however that purely frequential identification presents some well-known shortcomings. In particular, the so-identified symbol is in general ill adapted to the construction of efficient timerealizations. This is partly due to excessive numerical cost of quadrature approximations resulting from the intrinsic convolutive nature of the associated operator, sometimes with long memory (P. [START_REF] Bidan | Transient voltage distribution in motor windings fed by inverter: experimental study and modeling[END_REF] or even delay-like behaviors [START_REF] Montseny | Diffusive representation for operators involving delays, "applications of time-delay systems[END_REF]. Another heavy shortcoming is that frequential methods are incompatible with real-time identification (and so with pursuit when the symbol has the ability to evolve slowly). In opposite, time domain identification techniques do not present such drawbacks. However, their scope is in general not so wide. See for example (A. [START_REF] Monin | Arma lattice identification: a new hereditary algorithm[END_REF] for an interesting optimal method based on ARMA lattices.

Nevertheless, the notion of symbol can remain suitable for temporal identification, as it will be shown in this study. The proposed identification method is based on the so-called diffusive representation theory [START_REF] Montseny | Représentation Diffusive[END_REF] and was previously implemented, in a simpler version, on various problems (P. [START_REF] Bidan | Transient voltage distribution in motor windings fed by inverter: experimental study and modeling[END_REF][START_REF] Garcia | Identification of the dynamics of a lead acid battery[END_REF]. According to the diffusive representation approach, state formulations of diffusive type are available for a wide class of integral causal operators (including both rational and non rational ones) and the identification problem is carried over the so-called γ-symbol, defined as a suitable function (of an auxiliary real variable) easily deduced from the classical symbol. This γ-symbol entirely characterizes the associated operator, up to a suitable infinitedimensional state equation to be chosen a priori. In any case, cheap and precise finitedimensional approximations of this state representation can be straightforwardly built, running into finitedimensional problems easy to solve by means of standard techniques. Among the numerous advantages of this approach, we can mention in particular:

• a stable differential input-output time-formulation is available as soon as the γ-symbol is well-defined; • recursive identification algorithms are easily built and implemented, allowing real-time identification or even pursuit (in that case, the γ-symbol is depending on t); • similarly to purely frequential methods, no qualitative difference is made between rational and non rational components which can then be identified by the same process.

The paper is organized as follows. In section 2, we first introduce a few fundamental notions of the diffusive representation theory. In section 3, we describe the identification method under consideration, based on the previously introduced notions.

We then consider in section 4 a typical problem whose difficulty lies in the fact that both diffusive and delay components are involved but are no separable, and to which the diffusive representation approach is well-adapted: the identification of the acoustic impedance of an absorbent material devoted to aircraft motors noise reduction. Finally, a significant numerical example is given in section 5.

DIFFUSIVE FORMULATION OF CAUSAL INTEGRAL OPERATORS

In this section, we present a simplified version of a general methodology introduced and developed in [START_REF] Montseny | Représentation Diffusive[END_REF] in a general framework.

We consider a causal operator defined, on any continuous function u : R + → R, by

u → t 0 h(t -s) u(s) ds. ( 1
)
We denote H the Laplace transform of h and H(∂ t ) the convolution operator defined by (1).

Let u t (s) = 1 ]-∞,t] (s) u(s) the restriction of u to its past and u t (s) = u t (ts) the history of u.

From causality of H(∂ t ), we deduce:

H(∂ t )(u -u t ) (t) = 0 for all t;
then, we have for any continuous function u:

(H(∂ t ) u)(t)= L -1 (H Lu) (t)= L -1 H Lu t (t). ( 2
)
We define:

Ψ u (t, p) := e p t Lu t (p) = (Lu t ) (-p); (3)
by computing ∂ t Lu t , Laplace inversion and use of ( 2), it can be shown:

Lemma 1. 1. The function Ψ u is solution of the differential equation:

∂ t Ψ(t, p) = p Ψ(t, p) + u, t > 0, Ψ(0, p) = 0. 2. For any b 0, (H(∂ t ) u) (t) = 1 2iπ b+i∞ b-i∞ H(p) Ψ u (t, p) dp. (4)
We denote Ω the holomorphic domain of H (after analytic continuation). Let γ a closed1 simple arc in C -; we denote Ω + γ the exterior domain defined by γ, and Ω - γ the complementary of Ω + γ . By use of standard techniques (Cauchy theorem, Jordan lemma (M. [START_REF] Lavrentiev | Méthodes de la théorie des fonctions d'une variable complexe[END_REF]), it can be shown:

Lemma 2. For γ ⊂ Ω such that H is holomorphic in Ω + γ , if H(p) → 0 when p → ∞ in Ω + γ , then: (H(∂ t ) u) (t) = 1 2iπ γ H(p) Ψ u (t, p) dp, ( 5
)
where γ is any closed simple arc in Ω + γ such that γ ⊂ Ω - γ .

We now suppose that γ, γ are defined by functions of the Sobolev space2 W 1,∞ loc (R; C), also denoted γ, γ. We use the convenient notation μ, ψ = μ ψ dξ.

Under hypothesis of lemma 2, we have [START_REF] Montseny | Représentation Diffusive[END_REF]:

Theorem 3. If the possible singularities of H on γ are simple poles or branching points such that |H • γ| is locally integrable in their neighborhood, then:

1. with μ = γ 2iπ H • γ and ψ(t, .) = Ψ u (t, .) • γ: (H(∂ t ) u) (t) = μ, ψ(t, .) ; ( 6 ) 2. with 3 γn → γ in W 1,∞ loc and μ = γ 2iπ lim H • γn in the sense of measures: (H(∂ t ) u) (t) = μ, ψ(t, .) , ( 7
)
where ψ(t, ξ) is solution of the following evolution problem on (t, ξ) ∈ R * + ×R:

∂ t ψ(t, ξ) = γ(ξ) ψ(t, ξ) + u(t), ψ(0, ξ) = 0. (8) Definition 4. The measure μ defined in theo- rem (3) is called the γ-symbol of operator H(∂ t ).
Such a formulation can be extended to operators of the form

∂ t • H(∂ t )
where H(∂ t ) verifies the hypothesis of theorem (3). We indeed have (formally):

∂ t ψ(t, ξ) = γ(ξ) ψ(t, ξ) + u(t), ψ(0, ξ) = 0 ∂ t (H(∂ t ) u) (t) = μ, γ ψ(t, .) + μ, 1 u(t). ( 9
)
The state equation ( 8) is infinitedimensional. To get numerical approximations, we consider a discretization (ξ k ) k=1,n of the variable ξ and approximations μ n of the γ-symbol μ defined by atomic measures such as:

μ n = n k=1 μ k n δ ξ k . ( 10
)
Let us denote M n the space of atomic measures on the mesh {ξ k } k=1..n . If ∪ n M n is dense in the space of measures, we have [START_REF] Montseny | Représentation Diffusive[END_REF]:

μ n , ψ -→ n→+∞ μ, ψ ∀ψ.
So we have the finitedimensional approximate state formulation of

∂ t • H(∂ t ) : ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ ∂ t ψ(t, ξ k ) = γ(ξ k ) ψ(t, ξ k ) + u(t), k = 1 : n ψ(0, ξ k ) = 0, ∂ t (H(∂ t ) u) (t) n k=1 γ(ξ k ) μ k n ψ(t, ξ k ) + ( n k=1 μ k n ) u(t).
Note that when operator H(∂ t ) is not pseudodifferential, namely when some delay-type behaviors are present, the function γ necessarily verifies [START_REF] Montseny | Diffusive representation for operators involving delays, "applications of time-delay systems[END_REF]:

lim ξ→±∞ Im(γ(ξ)) = ±∞ and ∃a < 0, Re(γ(ξ)) a.
3 This convergence mode means that γn |K -γ |K → 0 and γ n |K -γ |K → 0 uniformly on any bounded set K.

OPTIMAL IDENTIFICATION OF THE γ-SYMBOL

In this section, we focus on the problem of identification of a convolution operator:

u -→ y = H(∂ t )u, (11)
where the symbol H(p) is supposed to verify the hypothesis of theorem 3. For simplicity, we present formal developments only; more details will be found in a further paper.

We denote by μ the γ-symbol of H(∂ t ), as defined in section 2. The problem then consists of identifying μ from the data u and y * , this last term denoting the measurement of y, in general sullied by some additive noise.

Given a suitable γ, we consider the diffusive realization (8,7) of H(∂ t ) and we denote ψ u the solution of (8). By defining A u the following linear operator:

A u : μ → μ, ψ u ,
we obviously have:

y = A u μ.
Given suitable Hilbert spaces E, F (not specified here), we consider the problem:

min μ∈E A u μ -y * 2 F ,
the solution of which is classically given by:

μ * = A † u y * (12)
where A † u denotes the pseudo-inverse of A u (A. [START_REF] Ben-Israel | Generalized inverses: theory and applications[END_REF]. So, in the sense of the hilbertian norm of F , the estimation μ * of μ is optimal.

From the numerical point of view, after suitable time discretization and approximation of μ such as described in section 2, we have:

(A u μ) i = n k=1 μ k n ψ(t i , ξ k );
then the operator A u can be expressed as a matrix (ψ u (t i , ξ k )) i=1:m,k=1:n and its pseudo-inverse is classically expressed4 (with m n):

A † u = (A * u A u + I) -1 A * u .
Remark 5. Recursive formulations of ( 12) can be established under the form (see (G. Garcia, 1998)):

μ * t = μ * t-Δt + K t-Δt (y * -A u μ * t-Δt ) |[0,t]
; such formulations allow real-time identification (or even the pursuit of μ in case of slowly varying operators Q(t, ∂ t )).

APPLICATION TO THE IMPEDANCE OF AN ABSORBENT MATERIAL

First, let us present the physical problem under consideration.

In order to reduce the noise generated by aircraft motors, the use of absorbent coverings is an efficient solution. For particular cases where stringent thermal or mechanical constraints are present, such as in hot zones of motors, porous metallic materials have been proposed and studied [START_REF] Gasser | Etude des propriétés acoustiques et mécaniques d'un matériau métallique poreux à base de sphères creuses de nickel[END_REF], (P.-A. [START_REF] Mazet | Control of aero-acoustic propagations with wall impedance boundary conditions: Application to a porous material model[END_REF]. Due to the complex structure of such materials, the generated acoustic impedances are non rational and complex functions of the frequency with unavoidable delay behaviors due to the presence of propagating waves in such media. In addition, their analytic determination is difficult and possibly poorly reliant. So, experimental measurements followed by identification processes remain necessary when accurate models are expected to carry out numerical simulations. In this way, it becomes possible, for a given material, to numerically compute the aeroacoustic field in a motor cavity equipped with absorbent walls in aim of evaluating the amount of noise reduction.

Fortunately, the impedance of an homogeneous 2D coating does not depend on the surface parameterization variables and therefore, at least for high frequencies, it can be identified from physical measurements relating to a small sample of the material (see figure 1). In the sequel, we show how the method introduced in section 4 can be implemented on such a problem. For simplicity of the analysis, we suppose that the problem described in figure 1 can be reduced to a 1D problem, which is legitimate if the source and the sensor are sufficiently far from the absorbent sample. This will be sufficient to introduce the principle of the method. For a more realistic analysis, the 3D problem can be similarly treated up to technical adaptations.

Under these conditions, the absorbent material sample can be summarized by an impedance op-erator Q(∂ t ) at point x = 1. The 1D problem of acoustic waves propagation can then be written as 5 :

∂ 2 t P = ∂ 2 x P + v ⊗ δ 0 on x ∈] -∞, 1], ( 13
)
where P is the acoustic pressure and v⊗δ 0 denotes the source at point x = 0, with v(t) a known timefunction. We have:

∂ x P |x=1 = Q(∂ t ) ∂ t P |x=1 . ( 14
)
Finally the acoustic pressure is measured at a point x m : y = P (t, x m ). ( 15)

Then, the identification problem is to build an approximation of operator Q(∂ t ) from the data y and v (note that in practice, the measurement y can be perturbed by an additive white noise, not explicitly specified here).

Let us compute the expression of the operator P (., 0) → y, denoted by K(∂ t ). For that, we must compute the reflexion coefficient at x = 1 for any frequency ω. Let ϕ an harmonic solution of ( 13) with v = 0. It can be expressed as ϕ = ϕ i + ϕ r where ϕ i denotes an incident wave (propagating in the direction of increasing x) and ϕ r the associated reflected wave (propagating in the sense of decreasing x), respectively given by: t+x) ,

ϕ i (t, ω, x) = e iω(t-x) , ϕ r (t, ω, x) = k(ω) e iω(
where k(ω) is the reflexion coefficient at frequency ω, depending on Q(iω). At point x = 1, the coefficient k(ω) verifies the impedance relation ( 14) which leads after simple computations to:

k(ω) = 1 + Q(iω) 1 -Q(iω) e -2iω .
So we have:

ϕ(t, ω, x) = (e -iωx + 1 + Q(iω) 1 -Q(iω) e iω(x-2) ) e iωt ,
and then (from 15): -2) .

K(iω) = e -iωxm + 1 + Q(iω) 1 -Q(iω) e iω(xm
On the other hand, the solution of ( 13) is classically given by:

P (t, x) = 1 2 (∂ -1 t v)(t -x)
; by harmonic synthesis, the input-output relation v → y can therefore be made explicit in the time domain (τ a denotes the translation operator

f (t) → f (t -a)): y = 1 2 τ xm • ∂ -1 t v + 1 2 M (∂ t ) • τ 2-xm • ∂ -1 t v (16)
with:

M (iω) := 1 + Q(iω) 1 -Q(iω) . ( 17
)
By defining the new data and notations:

ỹ := y -1 2 (τ xm + τ 2-xm ) • ∂ -1 t v, ṽ := 1 2 τ 2-xm v, M (∂ t ) := M (∂ t ) -1,
we obtain the suitable relation:

ỹ = ∂ -1 t • M (∂ t )ṽ (18)
on which we can directly apply the method presented in section 3 to identify the associated γsymbol.

Let μ denote indifferently the so-identified γsymbol or the exact γ-symbol of

∂ -1 t • M (∂ t )
. By considering the following relation deduced from (17):

Q(iω) = M (iω)
M (iω) + 2 and from elementary computations, we then deduce the following input-output state representation of operator

w → Q(∂ t ) w: ⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ ∂ t ψ = γ ψ - 1 2 γμ, ψ μ, 1 + 1 + 1 μ, 1 + 1 w, ψ 0 = 0 Q(∂ t )w = 1 2 γμ, ψ μ, 1 + 1 + μ, 1 μ, 1 + 1 w, (19)
which is exact (infinitedimensional) if μ is the exact γ-symbol, and approximate (finitedimensional6 ) if μ is the identified γ-symbol, of the form (10). Remark 6. A more precise analysis can be performed by replacing (13) by the 3D problem: ∂ 2 t P = ΔP + v ⊗ δ 0 . By use of Green functions, similar (but more technical) computations can be achieved and ( 18) is obtained again while the expression of ỹ and ṽ involves more complex (but explicitly known) operators derived from the Green functions.

NUMERICAL EXAMPLE AND CONCLUSION

We consider the porous material studied in (P.-A. [START_REF] Mazet | Control of aero-acoustic propagations with wall impedance boundary conditions: Application to a porous material model[END_REF], [START_REF] Gasser | Etude des propriétés acoustiques et mécaniques d'un matériau métallique poreux à base de sphères creuses de nickel[END_REF], the impedance of which has been analytically computed in (C. [START_REF] Casenave | State formulations of a porous medium model for boundary control of aeroacoustic waves[END_REF]. This impedance is expressed:

Q(iω) = ϕ χ eff (iω) ρ eff (iω) tanh iω χ eff (iω) ρ eff (iω) ,
where ϕ = 3.33 and the functions ρ eff (iω) and χ eff (iω) are respectively the so-called effective density of Pride et al. and the effective compressibility of Lafarge. These functions are expressed [START_REF] Gasser | Etude des propriétés acoustiques et mécaniques d'un matériau métallique poreux à base de sphères creuses de nickel[END_REF]:

ρ eff (iω) = eρ (1 + a (1+b iω) 1 2 iω ) χ eff (iω) = eχ (1 -c iω iω+a (1+b iω) 1 2
), with

ρ = ρ 0 α ∞ , χ = 1 P 0 , a = 8μ ρ 0 Λ 2 , a = 8μ ρ 0 Λ 2 , b = 1 2a , b = 1 2a , 0 < c = γ -1 γ < 1,
where e, ρ 0 , P 0 , μ, γ, α ∞ , Λ, Λ are physical parameters with values Λ = Λ = 0.1 10 -3 m, ρ 0 = 1.2 kg.m -3 , P 0 = 10 5 Pa, μ = 1.8 10 -5 kg.m -1 .s -1 , γ = 1.4, α ∞ = 1.3, e = 5 10 -2 m.

From standard analysis, it can be shown that the analytic continuation of Q presents singularities which can be asymptotically localized near a vertical straight line (see figure 2): operator Q(∂ t ) presents some underlying behaviors of delay type, inherited from propagative modes inside the porous material and from which result complex magnitude and phase variations of Q, making the identification problem possibly sensitive or even ill-posed. In aim of highlighting the efficiency of the previously introduced method, we consider the problem of numerical identification of such an operator Q(∂ t ) under its state formulation (19) in the ideal situation described in section 4.

The measurement data y * = y + ε have been generated on t ∈ [0, T ] from accurate numerical simulation of ( 16), with v a gaussian white noise and ε an output gaussian white noise with signalto-noise ratio equal to 26 dB.

The numerical parameters are Δt = 10 -5 s, T = 0.05 s, ξ min = ξ 1 = 10 3 rad/s, ξ max = ξ n = 3 10 5 rad/s, ξ k+1 = r ξ k , k = 1 : n -1, n = 100, γ(ξ) = |ξ| cos α + i ξ sin α, α = 100

• .

The comparison between the exact 7 and identified frequency responses relating to the impedance Q(iω) is available in figures 3 and 4. As expected, the identification is accurate in the frequency band [10 3 rad/s, 3 10 5 rad/s] covered by the set {γ(ξ k )} k=1:n . We can note that identification remains correct on the whole audio-frequency band, namely at low frequencies thanks to the asymptotically rational nature of Q(iω) which behaves like iω when ω → 0.

In figure 5, we can see the graph of the identified γ-symbol associated to operator ∂ -1 t • M (∂ t ) and from which is realized the operator Q(∂ t ) via the input-output formulation (19).

Beyond these quantitative results, we can add that this identification method is robust with respect to the level of the perturbation noise ε, the time of measurements T and the dimension n of the identified model. This last essential property is a consequence of the fact that the poles γ(ξ k ) of the identified transfer function are imposed a priori by the ξ-discretization of a suitable diffusive input-output model (namely ( 9)) mathematically well-posed and able to realize a wide class of operators in which belongs the non rational operator under consideration. It follows, thanks to the properties of diffusive representation, that the identified γ-symbol μ converges (in a sense not precised here), when n → +∞, α → π + 2 , ξ min → 0 + , ξ max → +∞ and T → +∞, to the exact one, and so for the identified impedance Q(iω) in the sense of the Fréchet space L 2 loc (R ω ).

Fig. 1 .

 1 Fig. 1. Physical measurement for the identification of the impedance of an homogeneous absorbing material.

Fig. 2 .

 2 Fig. 2. Singularities of the symbol Q(p).

Fig. 3 .7

 3 Fig. 3. Identified (-) and exact (---) and frequency responses of operator Q(∂ t ).

Fig. 4 .

 4 Fig. 4. Relative identification error |Q-Q * Q |.

Fig. 5 .

 5 Fig. 5. Identified γ-symbol μ(ξ), ξ > 0.

Possibly at infinity

W 1,∞ loc (R; C) is the topological space of measurable functions f : R → C such that f, f ∈ L ∞ loc (that is f and f are locally bounded in the almost everywhere sense).

As usual, is a small positive parameter devoted to numerical conditioning.

In that case, f, g = k f k g k .


Maximilien Guérin
email: maximilien.guerin@univ-paris3.fr

Network of Auxiliary Verb Constructions in Wolof: A constructional account of sentence structure

HAL is a

Introduction

In her reference monograph, [START_REF] Goldberg | Constructions: A Construction Grammar Approach to Argument Structure[END_REF] investigates one of the most central issue of grammar, i.e. the relation between verbs and sentential constructions. However, as noticed by Kuzar (2012: 1-4), while argument structure constructions play an important role in the formation of verbal sentences, they need a complementary linearization mechanism: the sentence pattern (S-pattern). In this paper, I will show how Wolof Predicative constructions provide data that question the concept of S-pattern.

1. Construction Grammar: some definitions 1.1. Construction. "C is a CONSTRUCTION iff def C is a form-meaning pair <F i , S i > such that some aspect of F i or some aspect of S i is not strictly predictable from C's component parts or from other previously established constructions." (Goldberg 1995: 4) (Croft 2001: 18) 1.2. Inheritance Links. "Constructions are specified as to which other, more abstract constructions they inherit from" (Goldberg 1995: 73), i.e. each construction is an instance of a more schematic construction (Figure 2). Moreover, "multiple inheritance paths are allowed" (ibid.), i.e. a construction may inherit from more than one other construction (Figure 3). (Croft 2001: 26) 1.3. Construction Networks. "Constructions form a structured inventory of a speaker's knowledge of the convention of their language. This structured inventory is usually represented by construction grammarians in terms of a taxonomic network of constructions. Each construction constitutes a node in the taxonomic network of constructions" (Croft 2001: 25) (Figure 4). However, because of multiple inheritance paths, grammatical constructions do not form a strict taxonomic hierarchy. Thus, I use the term 'extended network' to refer to a taxonomic network and the other constructions which dominate the constructions of that taxonomic network (Figure 5).

Linear order & Sentence patterns

Linear order is not specified in many of the more schematic constructions, even in a rigid word order language as English. For instance in ditransitive constructions, linear order is specified by other constructions (1a-d) (Croft 2001: 197).

( To give an account of linear order, Kuzar (2012) introduces a specific kind of sentential constructions (S-pattern), which constitutes a linearized list of component slots (2).

( (3) a. (i) They could easily find this book in the library.

) Verbal Sentence pattern in English (V S-pattern) NP [ VP V NP/PP NP/PP ] Sbj Pred Obj 1 Obj 2 I gave Becky a bottle ICCG8 -Maximilien Guérin -Network of Auxiliary Verb Constructions in Wolof 2/8 C1 C3 C2 C5 C7 C6 C8 C9 C10 C11 C4 C1 C3 C2 C5 C7 C6 C8 C9 C10 C11 C4 C12 C13 C14 C15 C16 2
V S-pattern, Declarative, Topic-Comment (ii) Could they easily find this book in the library?

V S-pattern, Interrogative, Topic-Comment (iii) This book -they could easily find it in the library.

V S-pattern, Declarative, Left-Dislocation (iv) This book -could they easily find it in the library?

V S-pattern, Interrogative, Left-Dislocation b. (i) It leaves me cold to see them dance together.

EV S-pattern, Declarative, Thetic (ii) Does it leave me cold to see them dance together?

EV S-pattern, Interrogative, Thetic (iii) What leaves me cold is to see them dance together.

EV sentence, Declarative, wh-cleft

Sentence Structure in Wolof

Wolof Language.

Wolof is a Niger-Congo language (belonging to the Atlantic family) mainly spoken in Senegal.

Verbal Predication.

In Wolof, verbal predication is based on a small number of constructions, the so-called "predicative constructions" (4a-h). These constructions are comparable to conjugations in that they express information about aspect (4a), tense (4b), mood (4c-d), or information structure (4e-h) [START_REF] Dialo | Structures verbales du wolof contemporain[END_REF][START_REF] Church | Le système verbal du wolof[END_REF][START_REF] Robert | Approche énonciative du système verbal : Le cas du Wolof[END_REF] These nine predicative constructions are in complementary distribution and form the basis of Wolof verbal predication, i.e. the speaker must choose one predicative construction to make a complete finite clause and the instantiation of a predicative construction matches a complete single clause.

With the exception of the imperative construction which is marked by an affix on the lexical verb (4d), Wolof predicative constructions display a mono-clausal split predicate structure involving a lexical verb and an auxiliary-like word, the predicative marker (PM). The PM carries the subject affix and the most part of TAM and information structure informations, and the lexical verb almost exclusively contributes lexical semantic content to the predication.

Linear Order.

In Wolof, S-patterns, TAM constructions and information structure constructions cannot be separated. Indeed, predicative constructions carry information dealing with tense, aspect, information structure or mood, and determine the linear order of all the sentential components (5-9).

(5) a. Perfect: IPFV-NEG-PST work 'You were not working.'

S V na-s o O b. Future: S dina-s o V O c. Optative: na-s o S V O d. Imperative: V-al-s o O e.

Sentence formation in Wolof.

Sentence formation in Wolof is based on the interaction of sentential constructions belonging to two major networks: the Network of argument structure constructions and the Network of predicative constructions (Figure 7). Some minor sentential constructions can mark the other dimensions: Negative, Past, Imperfective. Predicative constructions encode most constructional dimensions of the sentence, including linear order. Indeed sentence formation in Wolof is based on the interaction of sentential constructions belonging to two major networks: the Network of argument structure constructions and the Network of predicative constructions. Some minor sentential constructions can mark the other dimensions (polarity, imperfective aspect and past tense).

Figure 1 .

 1 Figure1. The symbolic structure of a construction(Croft 2001: 18)

Figure

  Figure 2. Instance Link

Figure

  Figure 4. Taxonomic Network Figure 5. Extended Network

Figure 6 .

 6 Figure 6. Argument structure mapped onto Sentence patterns

Figure 7 .

 7 Figure 7. Argument structure construction mapped onto AVC (Astu dóor na Aali 'Astou has hit Ali')

Figure 8 .

 8 Figure 8. Interaction of sentential constructions

  Just like Kuzar's S-patterns, every Wolof predicative construction has a default form, in which all the other dimensions are unmarked. By default, predicative constructions have a present perfective affirmative value (10a). To express imperfective aspect, the auxiliary verb di or its clitic form -y must be added to a predicative construction (10d). Negative polarity is generally encoded by the verb suffix -ul (10b), and past tense is usually encoded by the verb suffix -(w)oon (10c).

	Subject Focus: S/s-a o V O f. Presentative: S/s-angi o V O g. Compl. Focus: O/o la-s/S V h. Verb Focus: S da(fa)-s o V O S PM -s o V O a. Aali dina naan ndox Ali FUT drink water b. dina -a ko naan FUT 1SG 3SG drink S V PM -s o O a. Aali naan na ndox Ali drink PRF water b. naan na -a ko drink PRF 1SG 3SG O/o PM -s S V a. ndox la Aali naan water CFOC Ali drink b. lii la -a naan this CFOC 1SG drink PM -s o S V O a. na Aali naan ndox OPT Ali drink water b. na -a ko naan OPT 1SG 3SG drink da-ngeen liggéey PM-2PL work 'You have worked.' b. Negative da-ngeen liggéey -ul PM-2PL work -NEG 'You have not worked.' c. Past da-ngeen liggéey -oon PM-2PL work -PST 'You had worked.' d. Imperfective da-ngeen di liggéey 'You are working.' 3.4. Default form in Wolof. (10) a. Default form (6) (7) (8) (9) PM-2PL IPFV work	(Future) 'Ali will drink water.' 'I will drink it.' (Perfect) 'Ali has drunk water.' 'I have drink it.' (Complement focus) 'It is water that Ali has drunk.' 'It is this that I have drunk.' (Optative) 'May Ali drink water.' 'May I drink it.' e. Negative Past da-ngeen liggéey -ul-oon PM-2PL work -NEG-PST 'You had not worked.' f. Imperfective Negative da-ngeen d-ul liggéey PM-2PL IPFV-NEG work 'You are not working.' g. Imperfective Past da-ngeen d-oon liggéey PM-2PL IPFV-PST work 'You were working.' h. Imperfective Negative Past da-ngeen d-ul-oon liggéey PM-2PL

ICCG8 -Maximilien Guérin -Network of Auxiliary Verb Constructions in Wolof

Following the latest writing conventions, I write the presentative marker in two words: a ngi. Besides, the last vowel -i is a proximal deictic morpheme expressing a location close to the speaker. Two other morphemes may appear in this position: -a (distal) or -oog(u) (anaphoric)(Robert

2006: 164-165).ICCG8 -Maximilien Guérin -Network of Auxiliary Verb Constructions in Wolof

Abbreviations


research documents, whether they are published or not. The documents may come

On the design of observer-based fuzzy adaptive controller for nonlinear systems with unknown control gain sign


X Zhuang

M Lam

Chok Sing

C Dolabdjian

Investigation of the near-carrier noise for strain-driven ME laminates by using cross-correlation techniques

Keywords: Magnetoelectric, Equivalent magnetic noise, Strain modulation, Noise coherence

The near-carrier noise around the longitudinal mechanical resonance of a magnetoelectric laminated composite has been investigated. By simultaneously applying a high frequency electric field across the piezoelectric phase, the sensor response to low frequency magnetic signals can be shifted around the "carrier" frequency as side band modulation signals. This magnetoelectric response can appear either as an electric charge via piezoelectric to piezoelectric (PP) modulation effects or as a magnetic signal via piezoelectric to magnetic (PM) modulation effects. These two signals are detected either with a charge preamplifier or with a coil surrounding the sample and the low-frequency sensor response to the applied magnetic field can be recovered by using two independent synchronous detectors. We have designed an experimental setup to observe the direct (passive) low frequency noise and the noise corresponding to the two above modulations. Noise cross-correlating measurements were also carried out to investigate the origin of the near-carrier noise. No noise coherence was found between the direct low frequency noise and the noise resulting from either the PP or the PM modulations. However, a noise coherence factor of more than 50% has been found between the signals recovered from the two modulation techniques. A simple model has been used to explain this effect. The magnetoelectric sensor is considered as a nonlinear forced vibration system. Noise sources passing through such a system can be amplified and distributed around the carriers as side band noise where it hampers the equivalent magnetic noise performance. Electronic-thermal noise caused by dielectric dissipations in the piezoelectric phase can be considered as a noise source with a negligible contribution to the total noise floor. Mechanical-thermal low frequency excess noise is found to be the only intrinsic noise source which is filtered by the nonlinear ME system and is present as an output near-carrier noise which dominates the noise level after the demodulation processes.

I. INTRODUCTION

agnetoelectric (ME) sensors based on piezoelectric and magnetic thin films have huge potential for magnetic field detection via the product properties between magnetic and electric phases [START_REF] Fiebig | Revival of the magnetoelectric effect[END_REF], [START_REF] Nan | Multiferroic magnetoelectric composites: Historical perspective, status, and future directions[END_REF]. Magnetization and electric polarizations are mechanically coupled by mediated strains following certain direction. Magnetic sensibilities (usually defined as ME coefficients) are dominated by several aspects such as the piezoelectric and magnetic constants, a DC bias magnetic field, working frequencies and the volume fraction of the different phases [START_REF] Lam | Frequency response of magnetoelectric 1-3-type composites[END_REF], [START_REF] Zhou | Dielectric, magnetic, and magnetoelectric properties of laminated PbZr0.52Ti0.48O3/CoFe2O4/CoFe2 O4[END_REF]. The performance of a ME sensor can be enhanced by increasing its sensitivity. However, noise sources determine the limit of magnetic signal detection. Among several intrinsic noise sources, the one originating from the piezoelectric phase is predominant at frequencies of a few Hz [START_REF] Levinzon | Noise of Piezoelectric Accelerometer With Intergral FET Amplifier[END_REF]. In general, at low frequencies, electro-thermal noise associated to the delay effect between electric induction and electric field exceeds the contribution of the other noise sources in dielectric materials with a 1/f dependence of its voltage or current noise spectral density. In other words, this implies a harmful limit to the performance of a ME sensor at low frequencies and this is a physical obstacle in sensor designs and fabrication with piezoelectric dielectric materials, particularly with ferroelectric relaxors. Indeed, efforts on the reduction of such an electro-thermal noise can be achieved via advances of material science. But the progress of new materials cannot follow the rate of application development of such devices. However, another technique based on frequency shifting has been used to avoid the piezoelectric noise fluctuations at low frequencies. With the help of modulation techniques, low frequency magnetic signals can be transferred to higher frequencies around an exciting carrier. In this paper, the noise origins and distributions of such a technique are investigated and discussed via noise coherence measurements.

It has been demonstrated that when a modulation scheme is used, ME sensors have sensitivities which are different from those obtained under a passive sensing method [START_REF] Zhuang | Evaluation of Applied Axial Field Modulation Technique on ME Sensor Input Equivalent Magnetic Noise Rejection[END_REF] - [START_REF] Petrie | Enhancing the sensitivity of magnetoelectric sensors by increasing the operating frequency[END_REF]. In the passive mode, ME effects depend on the first order nonlinearity of the sensor and furthermore, the ME coefficient can be enhanced by applying an appropriate magnetic field bias. However, dominant perturbations originate from piezoelectric dissipations. Contrariwise, ME effects due to the second order nonlinearity are proportional to the magnitude of carriers when working under a modulating mode. This does not mean an infinite performance because of the limit from the electronics, the saturation of the sensor, etc... Besides, the observable noise for the modulation technique corresponds to the near carrier noise level which increases as the carrier amplitude is increased.

For the purpose of optimizing the noise reduction of ME sensors by such modulation techniques, the origin of near carrier noise must be investigated. If a ME sensor is excited by a harmonic electric carrier signal at high frequencies (e.g. several kHz) then, low frequency small magnetic signals (several Hz) can be detected across electrodes on piezoelectric layer. Thus, the PP modulation can be achieved. Unlike the PP modulation, another PM modulation can detect and recover the low frequency small magnetic signal via a pick-up coil surrounding the ME sensor.

In general, both the magnetic and electric modulation schemes can generate a strong strain in the ME sensors along the in-plane or along the out-of-plane directions. Thus, by applying an exciting electric carrier and measuring the noise M correlation from the piezoelectric and magnetic phases, the noise origins can be investigated.

II. EXPERIMENTAL DETAILS

Our sensors are made of commercial Metglas foils and piezoelectric fibers in order to form a 2-1 structure ME sample which has been demonstrated to produce a giant ME effect. Metglas foils of 25 µm in thickness were cut to a dimension of 8 cm long and 1 cm wide. Piezoelectric fibers were attached on both top and bottom surfaces by polymer isolators embedded with inter-digital electrodes by using epoxy resin to form a multi-push-pull structure [START_REF] Dong | Near-dieal magnetoelectricity in high-permeability magnétostrictives/piezofiber laminates with a (2-1) connectivity[END_REF].

Figure 1 shows our experimental setup. A magnetic flux density B inj generated by a pair of Helmholtz coils is used as a reference magnetic signal. A charge amplifier circuit with a 10 pF capacitor and a feedback resistance of 10 GΩ serves as a first stage amplifier electronics with a transfer gain of 10 10 V/C, as described in [START_REF] Zhuang | Evaluation of Applied Axial Field Modulation Technique on ME Sensor Input Equivalent Magnetic Noise Rejection[END_REF]. Generated charges can be recovered and transferred to this circuit via the amplifier inverting input and a voltage exciting the piezoelectric layers at 27.1 kHz can be applied simultaneously to the non-inverting input of charge amplifier via a generator V ex . The amplitude of this excitation signal has been chosen to 0.1 V to achieve an appropriate modulation effect without saturating the charge amplifier.

Thus, the reference signals under both PP and PM techniques can be recovered and amplified simultaneously via two independent electronic paths. The sensor response can be either observed directly in the passive (or classical) mode by connecting the output of the charge amplifier to Fil. 0, or detected after demodulation via Fil. 1. Fil. 0 is a low-pass filter which preserves the low frequency response to the reference signal while suppressing the high frequency components resulting from the excitation signal. These high frequency components are so high that they make the direct measurement impossible for the passive mode. Fil. 0 has a voltage gain of 100 and a bandwidth ranging from DC to 1 kHz. A simple 1 % attenuator was used to prevent the saturation of the filter. In order to recover the response of the ME sensor magnetic layers to the reference signal, a homemade 3000 turns pick-up coil surrounding the ME sensor was used, connected to a voltage amplifier having a voltage gain of 33 V/V. The output signal of this amplifier is then processed with Fil. 2 and then demodulated.

The band-pass filter Fil. 1 is used to select the components around the excitation frequency and it is associated with a compensation circuit, labeled "compensator" in figure 1, to reduce the large excitation carrier component, which is present at the output of the amplifier, to an appropriate level for the following demodulation stage. Filter Fil. 2 is identical to Fil. 1 for similar reasons. These two filters have the same gain of 100 and the same bandwidth from 3 kHz to 300 kHz to enhance the signals around the exciting carrier by filtering signals and perturbations at others frequencies far from the carrier. The noise and sensor response at low frequencies were then recovered via classical demodulation techniques by two phase-locked detectors.

III. THEORETICAL ANALYSIS

The limit of detection for a laminated magnetoelectric sensor at low frequencies, of the order of a few Hz, is often influenced by the electronic-thermal noise which is caused by the dielectric dissipation in the piezoelectric phase. Besides, another important intrinsic dissipation existing in the ME sensor is the mechanical-thermal noise. This noise originates from the molecule agitation and represents a white mechanical force noise at high frequencies [START_REF] Gabrielson | Mechanical-Thermal Noise in Micromachined Acoustic and Vibration Sensors[END_REF] and a 1/f noise dependence at quasi-static frequencies [START_REF] Riviere | Equation de Langevin et processus stochastiques[END_REF]. In the passive mode and at low frequencies (e.g. around 1 Hz), this mechanical noise is generally dominated by the electronicthermal noise in the piezoelectric material for a magnitude of several orders. However, the mechanical-thermal noise is partly involved in the nonlinear system of ME sensor. Thus, a frequency dependant noise contribution originating from the mechanical-thermal noise source appears around the carrier via the ME sensor strain nonlinearity. The magnetic, electric and mechanical coupling coefficient can be written as a following constitutive equation

  1 m m p p c c S S AT A e H A e E     (1)
where c and 1 c are linear and first order strain nonlinear rigidities, S and T are the strain and stress in sensor, A is the side-surface of sensor, A m and A p are respectively the sidesurface of the magnetic and piezoelectric layers, e m and e p are the magnetostrictive and piezoelectric constants, H and E are the magnetic and electric fields, respectively. The low frequency mechanical-thermal noise source results in a nearcarrier-noise around the exciting carrier. Through a magnetic and an electric detection, the noise formula can be written as H n (f) and E n (f) in [START_REF] Nan | Multiferroic magnetoelectric composites: Historical perspective, status, and future directions[END_REF]. It yields to

                1 0 1 0 2 2 p n n n c n c m m p m m n n n c n c p p p p e c E c H f S f S f f S f f e A c e A cE c E f S f S f f S f f e A c A                          (2)
where S n (f) is the mechanical noise source, c p the linear piezoelectric rigidity, E 0 the amplitude of the electric excitation signal, S n (f c -f) + S n (f-f c ) is the noise distribution near the excitation signal.

Several conclusions can be drawn: I) Any environmental magnetic noise can be amplified and distributed via the ME sensor nonlinearity as near-carrier noise. But, for our measurement setup, they are negligible in our magnetic shielded chamber. II) Electric-thermal dissipations are exhibited as an output noise which is not influenced by the nonlinear properties of the sensor. III) However, mechanicalthermal noise is partly influenced by the strain nonlinearity of the sensor and can be transferred to the output of the sensor as a phase noise distribution around the excitation carrier.

IV. RESULTS AND DISCUSSION

We first measured the sensor response to the reference magnetic signal at 2 Hz and the output noise floor in order to study later the increase of noise near the excitation carrier. We then applied a sinusoidal excitation voltage. The latter is modulated at the sensor output by the low frequency reference signal (cf. Fig. 2). We can obtain a giant enhancement in the sensitivity of ME effects by choosing the carrier frequency close to the mechanical resonance. The sensor response to the reference signal can be observed as two side-band signals at 27.1 kHz  2 Hz in Fig. 2. Signal spectra recorded for the PM modulation method are one order of magnitude higher than those given by the PP modulation method.

It is known that the sensitivity of a laminated ME sensor under nonlinear modulation techniques depend on several parameters namely the total nonlinearity, the piezoelectric and magnetic constants and the coupling coefficients. To better understand and compare the input equivalent magnetic noise floor corresponding to the classic passive mode, to the PM and PP modulation techniques, we performed the measurement of ME sensitivities of these three methods as a function of the frequency. We found that all of the three present good responses at low frequencies. As we mentioned above, the sensitivity of the PM modulation is higher than that for the PP modulation due to a higher nonlinear response. However, the low sensitivity in the classical passive mode can be explained by the absence of DC magnetic bias. In order to compare the noise level under the same magnetic conditions, no magnetic bias has been applied in all of our measurements.

Next, we obtained the transfer functions for the case of the passive mode, the PP and PM modulations. A reference magnetic flux of 38 nT rms was generated by a pair of commercial Helmholtz coils with a magnetic coefficient of 537 µT/A. The ME sensor, surrounded by a home-made pickup coil, is placed in the center of the Helmholtz coils by adjusting the sensing direction along the direction of the reference magnetic signal. These three transfer functions were directly measured at the V 0 , V 1 and V 2 outputs by using the "Sweep-sine" method of the spectrum analyzer.

Next step, we measured the noise floor with the three methods. We performed the experimental manipulations simultaneously two by two in order to analyze the coherence of the near-carrier noise. The sensor is driven by a voltage of 0.1 V in peak value around the sensor resonant frequency of f s = 27.1 kHz. By applying a reference magnetic signal to make sure that the system is operating correctly, the noise and the ME voltage introduced across inter-digital electrodes on both top and bottom sides of the piezoelectric layer of the laminate was measured as a function of frequency, using a dynamic spectrum analyzer. We measured the noise and the coherence for the passive mode and the two modulation techniques from either V 0 and V 1 or V 0 and V 2 outputs.

Figure 3 shows the input equivalent magnetic noise for the passive detection mode and the PP modulation technique and their coherence curve. It was found that the coherence is almost zero in the measurement bandwidth although the two noise curves were quite similar indicating that the noise sources are of different origins. Similar measurements have been performed between the passive mode and PM modulation. Because of the strong ME sensitivity related to the giant nonlinearity of magnetostrictive layers, the equivalent magnetic noise floor for the PM modulation technique is much lower as shown in Fig. 4. As for the preceding measurement, the coherence curve revealed different noise origins between the passive detection method and the PM modulation.

Finally, we determined the noise correlation between the PP and PM modulation techniques. We found a significant correlation factor (> 50%) at low frequencies as shown in Fig. 5.

These measurements results can be explained as follows: At low frequencies, the dominant noise for the passive mode originates from the electronic-thermal noise sources of sensor or noise sources in amplifier. None of them are real input noise sources for the ME nonlinear system which means that they cannot be distributed by the strong exciting carrier signal as a near-carrier noise. However, any environmental magnetic noise can be detected by the ME sensor system as an input perturbation. Such a noise source is detected by the ME sensor in the same way as the reference magnetic signal. In our magnetic shielded room, the environmental magnetic noise level is below the pT/√Hz at low frequencies. Therefore, the contribution from the environmental magnetic noise can be neglected compared to the observed noise level of the order of a few hundreds pT/√Hz at 1 Hz. Thus, the mechanical-thermal noise is the only noise source expected at low frequencies, modulated via the strain nonlinearity of laminated sensor.

V. CONCLUSIONS

Due to the intrinsic electric-thermal noise limit, caused by the dissipation in piezoelectric phases, the performance of ME sensors is hampered from achieving a considerable ME response under passive detection mode. Modulations techniques can be applied to overcome this inconvenience by shifting low frequency reference magnetic signals near a carrier at much higher frequencies via the nonlinearity of ME sensors. The intrinsic electric-thermal dissipation which presents itself as a dominant noise distribution in ME sensor can be avoided because of its presence after the nonlinear system of sensors. Therefore, the mechanical-thermal noise source in ME sensor appears as the only intrinsic noise source which is filtered by the nonlinear system and produces a frequency-dependent noise distribution around the exciting carriers. Our experiments with noise coherence techniques revealed the existence and the predominance of this noise distribution.

Fig. 2 .

 2 Fig. 2. Spectrogram of the near carrier noise and modulated signals around the excitation signal detected from the piezoelectric (red curve) and magnetic (blue curve) phases.

Fig. 1 .

 1 Fig. 1. Schematic setup of spectral cross-correlation technique for a ME sandwich laminates under the excitation of a sinusoidal electric polarization.

Fig. 3 .

 3 Fig. 3. Input equivalent magnetic noise for passive detection (black curve) and PP (red curve) modulation techniques and the spectral coherence between the two curves.

Fig. 4 .

 4 Fig. 4. Input equivalent magnetic noise for passive detection (upper black curve) and PM (lower blue curve) modulation technique and their spectral coherence (green curve).

Fig. 5 .

 5 Fig. 5. Input equivalent magnetic noise for PP (upper red curve) and PM (lower blue curve) modulation techniques and the spectral coherence (green curve) between these two curves.

ACKNOWLEDGEMENT

The authors acknowledge D. Viehland research group at Virginia Tech for providing the sensors.