Nicolas Chopin
email: nicolas.chopin@ensae.fr

Pierre Del Moral

Sylvain

Pierre Del Moral

S Rubenthaler

Stability of Feynman

published or not. The documents may come

Introduction

The most common application of the theory of Feynman-Kac formulae (see e.g. Del [START_REF] Del Moral | Feynman-Kac Formulae[END_REF]) is nonlinear filtering of a hidden Markov chain (Λ n ), based on observed process (Y n ). In such settings, the potential function at time n typically depends only on the current state Λ n . The uniform stability of the corresponding particle approximations can be obtained under appropriate conditions, see Section 7.4.3 of the aforementioned book and references therein. For a good overview of the theoretical and methodological aspects of particle approximation algorithms, also known as particle filtering algorithms, see also [START_REF] Doucet | Sequential Monte Carlo Methods in Practice[END_REF], [START_REF] Künsch | State space and hidden Markov models[END_REF], and [START_REF] Cappé | Inference in Hidden Markov Models[END_REF].

They are however several applications of practical interest where the potential function depends on the complete state trajectory Λ 0:n = (Λ 0 , . . . , Λ n ). The corresponding particle filtering algorithms still have a fixed computational cost per iteration, because the potential can be computed using recursive formulae. An important example is the class of conditional linear Gaussian dynamic models, where the conditioning is on some unobserved Markov chain Λ n . The corresponding particle algorithm is known as the mixture Kalman filter [START_REF] Chen | Mixture Kalman filters[END_REF], see also Example 7 in [START_REF] Doucet | On sequential Monte Carlo sampling methods for Bayesian filtering[END_REF]Andrieu and[START_REF] Andrieu | Particle filtering for partially observed Gaussian state space models[END_REF], for a related algorithm): the potential function at time n is then a Gaussian density, the parameters of which are computed recursively using the an inhomogeneous Markov chain, with initial probability distribution ζ, and transition kernel Q n . The observed process Y n admits Ψ n (y n |y 1:n-1 ; λ 0:n ) as a conditional probability density (with respect to an appropriate dominating measure) given Λ 0:n = λ 0:n and Y 1:n-1 = y 1:n-1 , where the short-hand v 0:n for any symbol v stands for the vector (v 0 , . . . , v n ). As explained in the Introduction, this quantity depends on the entire path λ 0:n , rather than the last state λ n . Following common practice, we drop dependencies on the y n 's in the notations, as the observed sequence y 0:n may be considered as fixed, and use the short-hand Ψ n (λ 0:n ) = Ψ n (y n |y 0:n-1 ; λ 0:n ). The model admits a Feynman-Kac representation which we describe fully in (2.1). We consider the following assumptions.

Hypothesis 1. For all n ≥ 1, the kernel Q n is mixing, i.e. there exists ε n ∈ (0, 1) such that

ε n ξ(A) ≤ Q n (λ n-1 , A) ≤ 1 ε n ξ(A)
for some ξ ∈ M + (E), and for any Borel set A ⊂ E, any λ n-1 ∈ E.

Hypothesis 2. For p large enough, and all n ≥ p, there exists a 'truncated' potential function Ψp n (λ n-p+1:n ) that depends on the last p states only, and that approximates Ψ n in the sense that

|Ψ n (λ 0:n ) -Ψp n (λ n-p+1:n )| ≤ φ n τ p Ψ n (λ 0:n ) ∧ Ψp n (λ n-p+1:n )
for some constants φ n and τ , φ n > 0, 0 < τ < 1, and all λ 0:n ∈ E n+1 . For convenience, we abuse notations and set Ψp n = Ψ n for p > n.

Hypothesis 3. There exists constants

a n , b n , n ≥ 0, a n ≥ 1, b n ≥ 1, such that 1 a n ≤ Ψ n (λ 0:n ) ≤ b n , 1 a n ≤ Ψp n (λ (n-p+1) + :n ) ≤ b n
for all λ 0:n ∈ E n+1 , using the short-hand k + = k ∨ 0 for any integer k.

The constants a n and φ n depend implicitly on the realisation y 1:n of the observed process. Hypotheses 1 and 3 are standard in the filtering literature; see e.g. Del [START_REF] Del Moral | Feynman-Kac Formulae[END_REF]. Hypothesis 2 formalises the fact that potential functions are computed using iterative formulae, and therefore should forget past states at an exponential rate. One may take Ψp n (λ n-p+1:n ) = Ψ n (z, . . . , z, λ n-p+1:n ) for instance, where z is an arbitrary element of E. We shall work out, in several models of interest, practical conditions under which Hypothesis 2 is fulfilled in Section 8.

We introduce the following notations for the forward kernels, for n ≥ 1:

γ n (λ 0:n-1 , dλ ′ 0:n ) = δ λ0:n-1 (dλ ′ 0:n-1 )Q n (λ n-1 , dλ ′ n )Ψ n (λ ′ 0:n )
where δ λ0:n-1 is the Dirac measure centred at λ 0:n-1 . The above kernels implicitly defines operators on measures and on test functions, i.e.,

γ n µ(f ) = γ n µ, f = µ(dλ 0:n-1 )γ n (λ 0:n-1 , dλ ′ 0:n )f (λ ′ 0:n ), for any µ ∈ M + (E n+1 ), any test function f : E n+1 → [0, 1], where M + (E k )
denotes the set of nonnegative measures w.r.t. E k , and P(E k ) the set of probability measures w.r.t. E k . We associate to γ n a "normalised" operator R n , such that, for any µ ∈ M + (E n ), R n µ is defined as:

R n µ(f ) = γ n µ(f ) γ n µ(1)
for any f : E n+1 → R + . Both the γ n 's and the R n 's may be iterated using the following short-hands, for 1 ≤ k ≤ n:

γ k:n µ = γ n . . . γ k µ, R k:n µ = R n . . . R k µ.
We have the following Feynman-Kac representation:

E(f (Λ 0:n )|Y 1:n = y 1:n ) = R 1:n ζ(f ) , (2.1)
∀n, ∀f : E n+1 → R + , where, as mentioned above, ζ the law of Λ 0 . Finally, we denote the total variation norm on nonnegative measures by • T V , the supremum norm on bounded functions by • ∞ , and the Hilbert metric by h(µ, µ ′ ) for any pair µ, µ ′ ∈ M + (E k ), k ≥ 1; see e.g. [START_REF] Atar | Exponential stability for nonlinear filtering[END_REF] or Le Gland and [START_REF] Le Gland | Stability and uniform approximation of nonlinear filters using the Hilbert metric and application to particle filters[END_REF], Definition 3.3. We recall that the Hilbert metric is scale invariant, and is related to the total variation norm in the following way, see e.g. Lemma 3.4 in Le Gland and [START_REF] Le Gland | Stability and uniform approximation of nonlinear filters using the Hilbert metric and application to particle filters[END_REF]:

µ -µ ′ T V ≤ 2 log 3 h(µ, µ ′ ) (2.2) h(Kµ, Kµ ′ ) ≤ 1 ε 2 µ -µ ′ T V (2.3)
provided K is a ε-mixing kernel. We can also derive the following properties from the definition of h (∀k

∈ N * , ∀µ, µ ′ ∈ M(E k )): ∀kernel Q, h(Qµ, Qµ ′ ) ≤ h(µ, µ ′ ) , (2.4) ∀nonnegative function ψ, h(ψµ, ψµ ′ ) ≤ h(µ, µ ′ ) (2.5)
with an equality in the latter equation if ψ is positive.

Local error induced by truncation

Until further notice, p is a fixed integer such that p ≥ 2 and such that Hypothesis 2 holds. Since our proofs involve a comparison between the true filter and a 'truncated' filter, we introduce the projection operator H p n which, for n ≥ p, associates to any measure µ(dλ 0:n ) ∈ M + (E n+1 ) its marginal w.r.t. its last p components, i.e. :

H p n (µ)(f ) = µ(dλ 0:n )f (λ n-p+1:n )
for any f : E p → R; for p > n, let H p n (µ) = µ. We also define the following 'truncated' forward kernels, for n ≥ p:

γp n (λ n-p:n-1 , dλ ′ n-p+1:n ) = δ λn-p+1:n-1 (dλ ′ n-p+1:n-1 )Q n (λ n-1 , dλ ′ n ) Ψp n (λ ′ n-p+1:n )
and the associated normalised operators, for µ ∈ M + (E p ), f :

E p → R + : Rp n µ(f ) = γp n µ(f ) γp n µ(1)
and set γp n = γ n , Rp n = R n for n < p. From now on, we will refer to the filter associated to these 'truncated' operators as the truncated filter.

We now evaluate the local error induced by the truncation.

Lemma 1. For all 1 ≤ k < n, and for all µ ∈ M + (E k ),

Rp k+1:n H p k R k µ -Rp k:n H p k-1 µ T V ≤ 2φ k τ p . Proof. Let f : E p∧(n+1) → [0, 1]. One has Rp k+1:n H p k R k µ(f ) = γp k+1:n H p k γ k µ(f ) γp k+1:n H p k γ k µ(1) Rp k:n H p k-1 µ(f ) = γp k:n H p k-1 µ(f ) γp k:n H p k-1 µ(1) where γp k+1:n H p k γ k µ(f ) = E n+1 µ(dλ 0:k-1 )Q k (λ k-1 , dλ k )Ψ k (λ 0:k )f (λ (n-p+1) + :n ) × n i=k+1 Q i (λ i-1 , dλ i ) Ψp i (λ (i-p+1) + :i ) and γp k:n H p k-1 µ(f ) = E n+1 µ(dλ 0:k-1 )Q k (λ k-1 , dλ k ) Ψp k (λ k-p+1:k )f (λ (n-p+1) + :n ) × n i=k+1 Q i (λ i-1 , dλ i ) Ψp i (λ (i-p+1) + :i ) hence γk+1:n H p k γ k µ(f ) -γk:n H p k-1 µ(f ) ≤ E n+1 µ(dλ 0:k-1 )Q k (λ k-1 , dλ k ) Ψ k (λ 0:k ) -Ψp k (λ (k-p+1) + :k ) ×f (λ (k-p+1) + :k ) n i=k+1 Q i (λ i-1 , dλ j ) Ψp i (λ (i-p+1) + :i ) ≤ φ k τ p E n+1 µ(dλ 0:k-1 )Q k (λ k-1 , dλ k )(Ψ k (λ 0:k ) ∧ Ψp k (λ (k-p+1) + :k )) ×f (λ (k-p+1) + :k ) n i=k+1 Q i (λ i-1 , dλ i ) Ψp i (λ (i-p+1) + :i ) ≤ φ k τ p γk+1:n H p k γ k µ(f ) ∧ γk:n H p k-1 γ k-1 µ(f )
f : f ∞ =1 E Rp k µ -Rk µ, f ≤ δ k
for some δ k ≥ 0, then, for all i ≥ 1 and µ ∈ P(E k∧p ),

sup f : f ∞=1 E Rp k:k+i µ -Rp k+1:k+i Rk µ, f ≤ 2(a k+1 . . . a k+i )(b k+1 . . . b k+i )δ k
where the expectation is with respect to the distribution of Rk .

Proof. Using the same ideas as above, one has, for f :

E (k+1-p+1)∧p → [0, 1], Rp k:k+i µ -Rp k+1:k+i Rk µ, f = γp k+1:k+i Rp k µ(f ) γp k+1:k+i Rp k µ(1) - γp k+1:k+i Rk µ(f ) γp k+1:k+i Rk µ(1)
.

In order to use inequality (3.6), compute

E( γp k+1:k+i Rp k µ(f ) -γp k+1:k+i Rk µ(f ) ) = E ( Rp k µ -Rk µ)(dλ (k-p+1) + :k ) k+i l=k+1 Q l (λ l-1 , dλ l ) Ψp l (λ (l-p+1) + :l )f (λ (k+i-p+1) + :k+i ) ≤ E b k+1 . . . b k+i ( Rp k µ -Rk µ)( f ) ≤ b k+1 . . . b k+i δ k where f is defined as f (λ (k-p+1) + :k ) = k+i l=k+1 Q l (λ l-1 , dλ l )f (λ (k+i-p+1) + :k+i ) ≤ 1.
and conclude by noting that γp k+1:k+i

Rp k µ(1) = ( Rp k µ)(dλ (k-p+1) + :k ) k+i l=k+1 Q l (λ l-1 , dλ l ) Ψp l (λ (l-p+1) + :l ) ≥ 1 a k+1 . . . a k+i since Rp k µ is a probability measure.

Mixing and contraction properties of the truncated filter

The truncated filter may be interpreted as a standard filter based on Markov chain Λp n = Λ (n-p+1) + :n . This insight allows us to establish the contraction properties of the truncated filter.

Lemma 3. One has:

h( Rp k+1:k+p µ, Rp k+1:k+p µ ′ ) ≤ 1 ε2 k+1,p µ -µ ′ T V and h( Rp k+1:k+p µ, Rp k+1:k+p µ ′ ) ≤ ρk+1,p h(µ, µ ′ ) where ε2 k,p = ε 2 k (a k . . . a k+p-2 )(b k . . . b k+p-2 ) , ρk,p = 1 -ε2 k,p 1 + ε2 k,p
, for all k ≥ 0, and all µ ,µ ′ ∈ P(E (k+1)∧p ).

Note εk,n must be interpreted as a mixing coefficient, and ρk,p as a Birkhoff contraction coefficient.

Proof. Using Hypothesis 3, one has:

Q k+p γk+1:k+p-1 µ = µ(dλ (k-p+1) + :k ) k+p i=k+1 Q i (λ i-1 , dλ i ) k+p-1 i=k+1 Ψp i (λ (i-p+1) + :i ) ≤ b k+1 . . . b k+p-1 µ(dλ (k-p+1) + :k ) k+p i=k+1 Q i (λ i-1 , dλ i ) ≤ b k+1 . . . b k+p-1 ε k+1 ξp (dλ k+1:k+p )
where ξp stands for the following reference measure:

ξp (dλ k+1:k+p ) = ξ(dλ k+1 ) k+p i=k+2 Q i (λ i-1 , dλ i ).
One shows similarly that

Q k+p γk+1:k+p-1 µ ≥ ε k+1 a k+1 . . . a k+p-1 ξp (dλ k+1:k+p ).
Hence kernel Q k+p γk+1:k+p-1 µ is mixing, with mixing coefficient εk+1,p . Following Lemma 3.4 in Le Gland and [START_REF] Le Gland | Stability and uniform approximation of nonlinear filters using the Hilbert metric and application to particle filters[END_REF],

h( Rp k+1:k+p µ, Rp k+1:k+p µ ′ ) = h(Q k+p γk+1:k+p-1 µ, Q k+p γk+1:k+p-1 µ ′ ) ≤ 1 ε2 k+1,p µ -µ ′ T V
using the scale invariance property of the Hilbert metric. Similarly, according to Lemma 3.9 in the same paper:

h( Rp k+1:k+p µ, Rp k+1:p µ ′ ) = h(Q k+p γk+1:k+p-1 µ, Q k+p γk+1:k+p-1 µ ′ ) ≤ 1 -ε2 k+1,p 1 + ε2 k+1,p h(µ, µ ′ ).

Propagation of truncation error

We establish first the two following lemmas.

Lemma 4. Let Rn : E n∧p → P(E (n+1)∧p ) be a sequence of (possibly random) probability kernels such that for all n ≥ 1 and µ ∈ P(E n∧p ), sup

f : f ∞=1 E Rp n µ -Rn µ, f ≤ δ n ,
where the expectation is w.r.t. the randomness of Rn , then, for all n ≥ 1 and all ζ ∈ P(E), one has

sup f : f ∞=1 E Rp 1:n ζ -R1:n ζ, f ≤ 8 log(3) n i=1   δ i ε2 i+1 ε2 i+p+1 ⌊ n-i p ⌋-1 j=2 ρi+jp+1,p  
where R1:n ζ = Rn . . . R1 ζ, and with the convention that empty products equal one.

Proof. The following difference can be decomposed into a telescopic sum:

Rp 1:n ζ -R1:n ζ = n i=1 Rp i+1:n Rp i R1:i-1 ζ -Rp i+1:n Ri R1:i-1 ζ .
We fix the integers i, n, and consider some arbitrary test function f . For i ≥ n -2p, one may apply Lemma 2:

sup f : f ∞ =1 E Rp i+1:n Rp i R1:i-1 ζ -Rp i+1:n Ri R1:i-1 ζ, f ≤ 2(a i+1 . . . a n )(b i+1 . . . b n )δ i ≤ 8 log(3) δ i ε2 i+1,p ε2 i+p+1,p since ε n ≤ 1, a n ≥ 1 and b n ≥ 1 for all n.
For i < n -2p, let k = ⌊(ni)/p⌋, then, using Lemma 3, Equations (2.2) to (2.5) one has 7) p. 160 of Le Gland and [START_REF] Le Gland | Stability and uniform approximation of nonlinear filters using the Hilbert metric and application to particle filters[END_REF], one gets

Rp i+1:n Rp i R1:i-1 ζ -Rp i+1:n Ri R1:i-1 ζ, f ≤ Rp i+1:n Rp i R1:i-1 ζ -Rp i+1:n Ri R1:i-1 ζ T V ≤ 2 log(3) h Rp i+1:i+kp Rp i R1:i-1 ζ, Rp i+1:i+kp Ri R1:i-1 ζ ≤ 2 log(3)ε 2 i+p+1,p × k-1 j=2 ρi+jp+1,p × Rp i+1:i+p ν -Rp i+1:i+p ν ′ T V where ν = Rp i R1:i-1 ζ, ν ′ = Ri R1:i-1 ζ. Applying (
Rp i+1:i+p ν -Rp i+1:i+p ν ′ T V ≤ 2 γp i+1:i+p ν -γp i+1:i+p ν ′ T V γp i+1:i+p ν(1)
.

where, using the same calculations as in Lemma 3,

γp i+1:i+p ν(1) ≥ ε i+1 a i+1 . . . a i+p and E γp i+1:i+p ν -γp i+1:i+p ν ′ T V = E x ′ ∈E p x∈E p (ν -ν ′ )(dx)γ p i+1:i+p (x, dx ′ ) ≤ sup x∈E p x ′ ∈E p γp i+1:i+p (x, dx ′ ) sup φ: φ ∞ =1 E(| ν -ν ′ , φ |) ≤ b i+1 . . . b i+p ε i+1 sup φ: φ ∞ =1 E(| ν -ν ′ , φ |)
which ends the proof.

Lemma 5. For all n ≥ 1 and all ζ ∈ P(E), one has

Rp 1:n ζ -H p n R 1:n ζ T V ≤ 4τ p log 3    n i=1 φ i ε2 i+1,p ⌊(n-i)/p⌋-1 j=1 ρi+jp+1,p   
with the convention that empty sums equal zero, and empty products equal one.

Proof. One has:

Rp 1:n ζ -H p n R 1:n ζ = n i=1 Rp i+1:n Rp i H p i-1 R 1:i-1 ζ -Rp i+1:n H p i R 1:i ζ For i ≤ n -p, let k = ⌊(n -i)/p⌋, then according to Lemma 3: Rp i+1:n Rp i H p i-1 R 1:i-1 ζ -Rp i+1:n Rp i+1 H p i R 1:i ζ T V ≤ 2 log 3 h Rp i+1:i+kp Rp i H p i-1 R 1:i ζ, Rp i+1:i+kp H p i R 1:i ζ ≤ 2 log(3)ε 2 i+1,p k-1 j=1 ρi+jp+1,p Rp i H p i-1 R 1:i-1 ζ -H p i R 1:i ζ T V
and ones concludes using Lemma 1. For i > np, one can apply Lemma 1 directly.

Coupling of particle approximations

We now introduce two interactive particle systems: the first particle system approximates the true filter, and is equivalent to the type of particle algorithms studied in this paper, and the second particle system approximates the truncated filter, and corresponds to an artificial algorithm that would not be implemented in practice. We work out a way of coupling both particle systems in order to evaluate the distance between the two (in a sense that is made clear below). We define, for n ≥ 0,

Qn,p λ (n-p) + :n-1 , dλ ′ (n-p+1) + :n = δ λ (n-p+1) + :n-1 (dλ ′ (n-p+1) + :n-1 ) ×Q n (λ n-1 , dλ ′ n ), Qn (λ 0:n-1 , dλ ′ 0:n ) = δ λ0:n-1 (dλ ′ 0:n-1 ) × Q n (λ n-1 , dλ ′ n ) . We define ∀ν ∈ M + (E n+1 ), ∀ measurable f : E n+1 → R + , ∀ν ′ ∈ M + (E p ), ∀ measurable g : E (n+1)∧p → R + , Ψ n .ν(f ) = ν(Ψ n f ) ν(Ψ n ) , Ψp n .ν ′ (g) = ν ′ ( Ψp n g) ν ′ (Ψ p n )
.

For any measurable space (E ′ , Ω ′ ) and any measure µ ′ ∈ P(E ′ ), we can take Z 1 , Z 2 , . . . i.i.d. of law µ ′ and define the random empirical measure, for N ≥ 1,

S N (µ ′ ) = 1 N N i=1 δ Zi .
Notice that, as the Z 1 , Z 2 , . . . are only given in law, we only define S N (µ) in law. We define the random operators R N n , Rp,N n (∀n) by: ∀µ ∈ P(E n ), R N n µ is a random weighted empirical measure such that

R N n µ = Ψ n .S N ( Qn µ) .
Similarly, ∀µ ′ ∈ P(E p∧n ), Rp,N n µ ′ is a random weighted empirical measure such that Rp,N n µ ′ = Ψp n .S N ( Qn,p µ ′ ) . (6.7)

As pointed above, R N n µ and Rp,N n µ ′ are only defined in law. Since ζ denotes the probability density of the first state Λ 0 , the particle system with N particles approximating the true filter at time n is defined by

R N n R N n-1 . . . R N 1 ζ
, and the particle system with N particles approximating the truncated filter at time n is defined by Rp,N n Rp,N n-1 . . . Rp,N 1 ζ. Lemma 6. There exists a coupling such that, for all k ≥ 1 and µ ∈ P(E k ):

sup f : f ∞≤1 E Rp,N k H p k-1 µ -H p k R N k µ, f ≤ φ k τ p .
As H p k R N k µ and Rp,N k H p k µ are defined to be random variables with such and such law, the term "coupling" means that we can define a random variable

(H k R N k µ, Rp,N k H k µ)
with the desired marginals. Proof. To prove the above result, we produce a coupling between the two random measures Rp,N k H p k-1 µ and

H p k R N k µ. Let Ψn (λ 0:n ) = Ψp n (λ (n-p+1) + :n ),
so that, for µ ∈ P(E k ), and using (6.7), one has

Rp,N k H p k-1 µ = H p k ( Ψk .(S N ( Qk µ))
in the sense that both sides define the same distribution. Let χ 1 , . . . , χ N i.i.d. ∼ µ Qk , where χ i is a vector λ 0:k,i , for i = 1, . . . , N , and χi denotes its projection on the p last components, χi = λ (k-p+1) + :k,i , then

1 N j=1 Ψ k (χ j ) N i=1 Ψ k (χ i )δ χi has same law as H p k R N k µ and 1 N j=1 Ψk (χ j ) N i=1
Ψk (χ i )δ χi has same law as Rp,N k H p k-1 µ .

For any f such that f ∞ ≤ 1 (using a classical result on empirical measures):

| Rp,N k H p k-1 µ -H p k R N k µ, f | ≤ 1 2 N i=1 Ψ k (χ i ) N j=1 Ψ k (χ j ) - Ψk (χ i ) N j=1 Ψk (χ j ) ≤ 1 2 N i=1   Ψ k (χ i ) -Ψk (χ i ) N j=1 Ψ k (χ j ) + Ψk (χ i ) N j=1 ( Ψk (χ j ) -Ψ k (χ j )) N j=1 Ψ k (χ j ) N j=1 Ψk (χ j )   ≤ φ k τ p 2 N i=1   Ψ k (χ i ) N j=1 Ψ k (χ i ) + Ψk (χ i ) N j=1 Ψk (χ i ) ∧ Ψ(χ i ) N j=1 Ψ(χ j ) N j=1 Ψ(χ j )   ≤ φ k τ p
using Hypothesis 2, from which we deduce the result.

Main result

We are now able to derive estimates of the error

E p n,N (y 1:n ) = sup f : f ∞=1 E N H p n R 1:n ζ -H p n R N 1:n ζ, f |Y 0:n = y 0:n
induced by the particle approximation of the true filter, for the marginal filtering distribution of the p last states, provided p ≤ n. The expectation E N is with respect to the randomness of the N particles, and the functions f are E (n+1)∧p → R. Note that E p n,N (y 1:n ) is by construction an increasing function of p.

Theorem 7. For any ζ ∈ P(E), and any test function w.r.t. E (n+1)∧p ,

E p n,N (y 1:n ) ≤ 4 log 3 n i=1 δ i ε2 i+1,p ε2 i+p+1,p ⌊(n-i)/p⌋-1 j=2 ρi+jp+1 (7.8)
where

δ i = 3τ p φ i + 4a i b i √ N .
Proof. We first study the following local error, for

µ ∈ P(E n ), sup f : f ∞ =1 E N Rp n H p n-1 µ -H p n R N n µ, f Y 0:n = y 0:n
where the difference of operators can de decomposed into:

Rp n H p n-1 -H p n R N n = Rp n H p n-1 -Rp,N n H p n-1 + Rp,N n H p n-1 -H p n R N n .
To bound the first term, one may use (25) p. 162 of Le Gland and [START_REF] Le Gland | Stability and uniform approximation of nonlinear filters using the Hilbert metric and application to particle filters[END_REF], for ν = H p n-1 µ and Hypothesis 3:

E N Rp n ν -Rp,N n H p n-1 ν, f ≤ 2a n b n √ N
and, for the second term, one may apply Lemma 6:

E N Rp,N n H p n-1 µ -H p n R N n µ, f ≤ φ n τ p so that sup f : f ∞ =1 E N Rp n µ -H p n R N n µ, f ≤ δ ′ n for δ ′ n = 2a n b n / √ N + φ n τ p .
This local error is propagated using Lemma 4:

E N Rp 1:n ζ -H p n R N 1:n ζ, f ≤ 8 log(3) n i=1   δ ′ i ε2 i+1 ε2 i+p+1 ⌊ n-i p ⌋-1 j=2 ρi+jp+1,p   .
To conclude, one may decompose the global error as follows:

H p n R 1:n ζ -H p n R N 1:n ζ = H p n R 1:n ζ -Rp 1:n ζ + Rp 1:n ζ -H p n R N 1:n ζ .
where the second term is bounded above, and the first term is directly bounded using Lemma 5.

Since p is an arbitrary parameter, one may minimise the error bound with respect to p. For instance, one has the following result for time-uniform estimates. As noted above, the error E p n,N (y 1:n ) is an increasing function of p, so the bound below applies a fortiori to E 1 n,N (y 1:n ), the particle error corresponding to the marginal filtering distribution of the last state Λ n .

Corollary 8. If there exists constants c, ε, φ > 0 such that, almost surely, a n b n ≤ c, ε n ≥ ε, and φ n ≤ φ, then, provided τ c 3 < 1, the particle error is bounded almost surely as follows:

E p n,N (y 1:n ) ≤ C {log(N ) + D} 1 √ N 1+3 log c/ log τ ,
for N large enough, where

C = 16 ε 6 c 2 -1 log τ 4c 3φ 3 log c/ log τ , D = 2 log(3φ/4cτ ), and
p = log 4c 3φ √ N / log τ . (7.9)
Proof. Under these conditions, the RHS of (7.8) is smaller than or equal to:

E p n,N (y 1:n ) ≤ 4 log 3 c 2(p-2) ε 4 3φτ p + 4c √ N n i=1 1 -ε 2 c -(p-2) ⌊(n-i)/p⌋-2 (7.10) ≤ 4 log 3 c 2(p-2) ε 4 3φτ p + 4c √ N n-1 i=0 1 -ε 2 c -(p-2) i/p-1 ≤ 4 log 3 c 2(p-2) ε 4 3φτ p + 4c √ N 1 -ε 2 c -(p-2) -1 1 -1 -ε 2 c -(p-2) 1/p ≤ 4c 3(p-2) ε 6 3φτ p + 4c √ N p
for p large enough, since (1x) a ≤ 1ax for a ∈ (0, 1), x ∈ (0, 1), so, provided c 3 τ < 1, one may take p as in (7.9), which gives:

E p n,N (y 1:n ) ≤ 32 ε 6 c 2 4c 3φ 3 log c log τ log N + 2 log(3φ/4c) -2 log τ + 1 1 √ N 1+ 3 log c log τ
and conclude.

Obviously, this is a qualitative result, in that there are many practical models where such time-uniform, deterministic bounds are not available. For specific models, one may be able instead to use (7.8) in order to establish the asymptotic stability of the expected particle error, where the expectation is with respect to observed process (Y n ). We provide an example of this approach in Section 8.

Applications to practical models

In this section, we apply our general result to three practical models. We keep the same settings and notations, i.e. the observed process (Y n ) admits some probability distribution conditional on the path Λ 0:n = λ 0:n of a Markov chain (Λ n ), with initial distribution ζ and Markov transition Q n , which fulfil Hypothesis 1, see Section 2. We derive conditions on the model parameters that ensure asymptotic stability of the particle error; in particular, these conditions imply that Hypotheses 2 and 3 are verified.

We state the following trivial result for further reference. Let (f, g) a pair of probability densities (f, g) on E, then:

∀x ∈ E, | log {f (x)} -log {g(x)} | ≤ c ⇒ ∀x ∈ E, |f (x) -g(x)| ≤ (e c -1) {f (x) ∧ g(x)} (8.11)
for c ≥ 0.

GARCH Mixture model

We assume that the observed process is such that

Y n = σ n (Λ 0:n )Z n , n ≥ 1,
where the Z n 's are i.i.d. N (0, 1) random variables, and the variance function σ 2 n is defined recursively, for n ≥ 1:

σ 2 n (λ 0:n ) = α(λ n ) + β(λ n )Y 2 n + γ(λ n )σ 2 n-1 (λ 0:n-1 ) (8.12)
and σ 2 0 (λ 0 ) = α(λ 0 )/ {1γ(λ 0 )} , where α, β and γ are E → R + functions. Conditional on Λ 0:n , (Y n ) is a GARCH (generalised autoregressive conditional heteroskedasticity) process [START_REF] Bollerslev | Generalized autoregressive conditional heteroskedasticity[END_REF]; see [START_REF] Chopin | Dynamic detection of change points in long time series[END_REF] for a finance application of such a GARCH mixture model.

The potential functions equal

Ψ n (λ 0:n ) = 1 2πσ 2 n (λ 0:n ) exp - y 2 n 2σ 2 n (λ 0:n ) ,
for λ 0:n ∈ E n+1 , and (Λ n ) is a Markov process, with Markov kernels Q n , which satisfy Hypothesis 1. The functions α, β and γ are assumed to be bounded as follows:

0 < α min ≤ α(λ) ≤ α max , 0 ≤ β min ≤ β(λ) ≤ β max < 1, 0 ≤ γ min ≤ γ(λ) ≤ γ max < 1.
We first consider the case where β(λ) = 0 for all λ ∈ E. As mentioned in the introduction, this simplified model can be interpreted as a standard hidden Markov model, with observed process (Y n ), and Markov chain (κ n ) = Λ n , σ 2 n (Λ 0:n ) . However, since σ 2 n (Λ 0:n ) is a deterministic function of σ 2 n-1 (Λ 0:n-1 ) and λ n , it does not have mixing or similar properties that are usually required to obtain estimates of the particle error. Instead, analysing this model as a Feynman-Kac flow with iterative, path-dependent potential functions make it possible to derive such estimates.

Lemma 9. For the simplified model described above (with β = 0), the expected particle error of the corresponding particle approximation is uniformly stable in time, i.e. there exists constants C, D, such that

E E p n,N (Y 1:n ) ≤ C {log(N ) + D} 1 √ N 1+3 log c/ log τ ,
where p is given by (7.9), provided ι < 2 and τ c 3 < 1, where τ = γ max , c = (2/ι -1) -1/2 , and

ι = α max (1 -γ min ) α min (1 -γ max ) .
Proof. From (8.12), one sees the process

σ 2 n is bounded, σ 2 min ≤ σ 2 n (λ 0:n ) ≤ σ 2 max
for all λ 0:n ∈ E n+1 , where

σ 2 min = α min 1 -γ min , σ 2 max = α max 1 -γ max .
so, for a given sequence observations y 1:n , Hypothesis 3 is verified with:

1 a n = 1 2πσ 2 max exp - y 2 n 2σ 2 min , b n = 1 2πσ 2 min exp - y 2 n 2σ 2 max ,
provided the truncated potential is taken as:

Ψp n (λ n-p+1:n ) = Ψ n (z, . . . , z, λ n-p+1:n )
where z is an arbitrary element of E. For Hypothesis 2, one has, for any

λ 0:n , λ ′ 0:n ∈ E (n+1) such that λ (n-p+1) + :n = λ ′ (n-p+1) + :n : |log Ψ n (λ 0:n ) -log Ψ n (λ ′ 0:n )| ≤ 1 2 log σ 2 n (λ 0:n ) -log σ 2 n (λ ′ 0:n ) + y 2 n 2 1 σ 2 n (λ 0:n ) - 1 σ 2 n (λ ′ 0:n ) ≤ σ 2 min + y 2 n 2σ 4 min σ 2 n (λ 0:n ) -σ 2 n (λ ′ 0:n )
where σ 2 n is contracting, in the sense that, for n ≥ p,

σ 2 n (λ 0:n ) -σ 2 n (λ ′ 0:n ) = p-1 i=0 γ(λ n-i ) σ 2 n-p (λ 0:n-p ) -σ 2 n-p (λ ′ 0:n-p ) ≤ 2γ p max σ 2 max .
Thus, using (8.11), and the fact that (e x -1)/x is an increasing function, Hypothesis 2 is verified with τ = γ max and

φ n = τ -q exp τ q σ 2 max σ 2 min + y 2 n σ 4 min -1 ,
for any q ≤ p. Finally, to compute the expectation with respect to process (Y n ) of the error bound (7.8), one may use repetitively the following results:

E exp aY 2 n |Y 1:n-1 ≤ 1 -2aσ 2 max -1/2
for a < 1/2σ 2 max , using standard calculations and the fact that Y n , conditional on Y 1:n-1 and Λ 0:n = λ 0:n is N 0, σ 2 n (λ 0:n ) . This implies in particular that:

E [a n b n |Y 1:n-1 ] ≤ 2 σ 2 min σ 2 max -1 -1/2 = c
where the constant c is well-defined since σ 2 max /σ 2 min < 2, then by Jensen inequality,

E 1 a n b n Y 1:n-1 ≥ c -1 ,
and similarly,

E [φ n |Y 1:n-1 ] ≤ τ -q exp τ q σ 2 max σ 2 min 1 -2τ q σ 4 max σ 4 min -1/2 -1 = φ
where φ is properly defined for q large enough. Using the above results recursively on the sum on the RHS of (7.8), one obtains the same expression as in (7.10) for the error bound than in Corollary 8 for time-uniform estimates (with the values of c, φ, τ as defined above), and concludes similarly.

If β is allowed to take positive values, stability results may be obtained under more restrictive conditions. In particular, one may impose that γ is a constant function.

Lemma 10. For the general mixture GARCH model defined above, the expected particle error is uniformly stable in time, i.e. there exist constants C, D, such that

E E p n,N (Y 1:n ) ≤ C {log(N ) + D} 1 √ N 1+3 log c/ log γ
provided γ is a constant function, γ(λ) = γ, τ c 3 < 1, ϑ < 2, where τ = γ, c = (2/ϑ -1) -1/2 , p is given by (7.9), and

ϑ = α max α min ∨ β max β min .
Proof. We follow the same lines as above, except that the bounds of the process σ 2 n (λ 0:n ) must be replaced by:

σ 2 min (n) = γ n 1 -γ α min + n-1 k=0 (α min + β min y 2 n-k )γ k , σ 2 max (n) = γ n 1 -γ α max + n-1 k=0 (α max + β max y 2 n-k )γ k ,
which, by construction, are such that

σ 2 max (n) σ 2 min (n) ≤ ϑ < 2.
Hence, one has again

E [a n b n |Y 1:n-1 ] ≤ 2 σ 2 min σ 2 max -1 -1/2 = c
and that, for λ 0:n , λ ′ 0:n such that λ n-p+1:n

= λ ′ n-p+1:n , |µ n (λ 0:n ) -µ n (λ ′ 0:n )| (8.19) ≤ M n-1   p i=1 hi a n-i i-1 j=1 (1 -a n-j ) -a ′ n-i i-1 j=1 (1 -a ′ n-j ) + 2 hp+1  
where a n-i , a ′ n-i are short-hands for a n-i (λ 0:n-i ), a n-i (λ ′ 0:n-i ). The sequence a n itself in contracting, since, from (8.15), one has, for i < p:

a n-i -a ′ n-i ≤ v σ 4 σ 2 n-i (λ 0:n-i ) -σ 2 n-i-1 (λ ′ 0:n-i ) ≤ vC σ σ 4 τ p-i σ so (8.19) and the fact that |xy -x ′ y ′ | ≤ |x -x ′ | + |y -y ′ | provided x, x ′ , y, y ′ ∈ [0, 1] leads to |µ n (λ 0:n ) -µ n (λ ′ 0:n )| ≤ M n-1 vC σ σ 4 p i=1 hi τ p-i σ + . . . + τ p-1 σ + 2 hp+1 ≤ M n-1 vC σ τ p-1 σ σ 4 p i=1 hi τ -i σ -1 τ -1 σ -1 + 2 hp+1 ≤ M n-1 C µ τ p
for τ = τ σ ∨ h, and a well chosen value of C µ .

We are now able to state the main result.

Lemma 13. For the model above, the particle error is bounded uniformly in time, i.e. there exist C, D, such that

E p n,N (y 1:n ) ≤ C {log(N ) + D} 1 √ N 1+3 log c/ log τ ,
almost surely, for p given by (7.9), provided the realizations y n are bounded, i.e. |y n | ≤ C y for all n ≥ 1, and that τ c 3 < 1, with τ = h ∨ τ σ and

c = σ σ exp C 2 y σ 2 1 + āh 1 -ãh 2 , τ σ = 1 1 + w/v + 2 w/v + w 2 /v 2 < 1.
Proof. This proposition is a direct application of Corrolary 8, so we need only to prove that Hypotheses 2 and 3 are fulfilled. For Hypothesis 2, one may take

1 a n = 1 √ 2πσ 2 exp - C 2 y σ 2 1 + āh 1 -ãh 2 , b n = 1 2πσ 2
so that a n b n ≤ c for c defined above. For Hypothesis 3, one has:

2 |log Ψ n (λ 0:n ) -log Ψ n (λ ′ 0:n )| ≤ log σ 2 n (λ 0:n ) -log σ 2 n (λ ′ 0:n ) + {y n -µ n (λ 0:n )} 2 σ 2 n (λ 0:n ) - {y n -µ n (λ ′ 0:n )} 2 σ 2 n (λ ′ 0:n )
Clearly, the reformulated model satisfies Hypothesis 1: the (Λ n ) are i.i.d., hence they form a Markov chain with mixing coefficient ε n = 1. If we assume that the Ψ n (λ 0:n ) are such that Hypotheses 2 and 3 hold as well, then we can apply directly Theorem 7. However, the path-dependent formulation of this model is artificial, and, in practice, we are interested in filtering the process X n , conditional on the Y n 's, rather than filtering the Λ n 's, again conditional on the Y n 's. More precisely, we wish to approximate the conditional expectation of

g(X n ) = g   p-1 k=0 ρ k λ n-k + n k=p ρ k λ n-k   ,
for some bounded function g, and, provided g is also Lipschitz, with constant K, and that the λ n 's lie in interval [-l, l], for some l ≥ 0, one has:

g   p-1 k=0 ρ k λ n-k + n k=p ρ k λ n-k   -g p-1 k=0 ρ k λ n-k ≤ Kl 1 -τ τ p ,
where τ = |ρ|. Therefore, we must consider an additional term in the particle error attached to the filtering of (X n ), which stems from the difference between the filtering distribution of X n and that of Λ n-p+1:n , for some integer p. Consider the following estimate of the particle error for functions of X n :

E X n,N (y 1:n ) = sup g: g ∞ =1,g∈Lip(K) E N R 1:n ζ -R N 1:n ζ, f g |Y 0:n = y 0:n
where Lip(K) denotes the set of Lipschitz functions with Lipschitz constant K, and f g is the function E n+1 → R such that

f g (λ 0:n ) = g n k=0 ρ k λ n-k ,
i.e., loosely speaking, f g (λ 0:n ) = g(x n ), where x n must be substituted by its expression as a function of λ 0:n .

Lemma 14. For the state-space model described above, one has, for any n ≥ p,

E X n,N (y 1:n ) ≤ E p n,N (y 1:n ) + Kl 1 -τ τ p .
Taking into account this additional error term, we can derive time-uniform estimates of the stability of the particle algorithm. For the sake of space, we focus on the following simple example: Y n ∈ {-1, 1}, Y n = 1 with probability 1/(1 + e Xn ), Y n = -1 otherwise. The potential function (for the model in its standard formulation) equals:

Ψ X n (x n ) = 1 1 + e ynxn .
We recall that the support of the (Λ n ) is [-l, l], and therefore X n ∈ [-l ′ , l ′ ] almost surely, with l ′ = l/(1-τ ). Thus, Hypothesis 3 holds for b n = 1/(1+e -l ′ ), The above model does not fulfil the usual conditions required in standard stability results, see e.g. Del Moral (2004, Section 7.4.3), because the Markov chain (X n ) is not mixing. Thus, it is remarkable that the time-uniform stability of this model is established using a Feynman-Kac formulation with path-dependent potentials.

Conclusion

To extend our results to a broader class of models, three directions may be worth investigating. First, it may be possible to bound directly the particle error, without resorting to a comparison with an artificial, truncated potential function. It seems difficult however to avoid some form of truncation, as the path process Λ 0:n itself does not benefit to any sort of mixing property, while fixed segments Λ n-p+1:n do. Second, one may try to loosen Hypothesis 1 (Markov kernel is mixing) and Hypothesis 3 (potential function is bounded), using for instance [START_REF] Oudjane | Stability and uniform particle approximation of nonlinear filters in case of non ergodic signals[END_REF]'s approach. Third, it seems possible to adapt our general result on the particle error bound to several models not considered in this paper, in particular standard models with potential functions depending on the last state only, by using and extending the approach developed in the previous Section.

  a n = 1 + e l ′ . For Hypothesis 2, standard calculations show that, for two vectors λ 0:n and λ ′ 0:n such that λ n-p+1:n = λ ′ n-p+1:n , one has|log Ψ n (λ 0:n )log Ψ n (λ ′ 0:n )| ≤ n k=p ρ k (λ n-kλ ′ n-k ) ≤ 2l ′ τ p provided τ = |ρ|.Hence, using (8.11) inequality, Hypothesis 2 holds, withφ n = e 2l ′ -1.For this specific model, we have the following result.Lemma 15. For the specific model described above, and provided cτ 3 < 1, where τ = |ρ|, c = e l ′ , one has:E X n,N (y 1:n ) ≤ C {log(N ) + D} 1 √ N 1+3 log c/ log τ + E √ Nwhere C and D were defined in Corollary 8, φ = e 2l ′ -1, and E = 4Kl ′ c/3φ.

Acknowledgements

Part of this work was performed while the authors were invited by the SAMSI institute during the 2008-2009 Program on Sequential Monte Carlo Methods.

and the rest of the calculation is identical to those of previous Lemma, with τ = γ.

Mixture Kalman model

We focus on an univariate linear Gaussian model, i.e. conditional on Markov process (Λ n ), one has X 0 = 0 almost surely, and, for n ≥ 1,

where the V n 's and the W n 's are independent N (0, 1) variables, and h, v, w are real-valued functions. Using the recursions of the Kalman-Bucy Filter [START_REF] Kalman | New results in linear filtering and prediction theory[END_REF], one is able to marginalise out the process X n , and compute recursively the probability density of Y n , conditional on Λ 0:n = λ 0:n , in the following way:

where, the following quantities are defined recursively: for n ≥ 1,

and m 0 (λ 0 ) = c 0 (λ 0 ) = 0. We make the following assumptions:

1. Functions v and w are bounded as follows: for all λ ∈ E,

2. Function h is bounded as follows: for all λ ∈ E,

We first prove the following intermediate results.

Lemma 11. The sequence σ 2 n is bounded and uniformly contracting, i.e. for all p ≥ 1, for all λ 0:n , λ ′ 0:n , such that λ n-p+1:n = λ ′ n-p+1:n , one has

where σ 2 = v + w, σ2 = ( h2 + 1)v + w, C σ = h2 v/τ σ , and

Proof. From (8.17), one deduces that

In addition, (8.18) implies that

.

It is easy to show that, for a fixed λ, the derivative of Υ (c, λ) with respect to c is bounded from above by τ σ as defined above. Thus, Υ (c, λ) is a contracting function, and, by induction, for n ≥ p,

where τ σ and C σ were defined above.

Lemma 12. The sequence µ n is bounded and contracting in the sense that there exists C µ > 0 such that, for all p ≥ 1, for all n ≥ p, and λ 0:n , λ ′ 0:n , such that λ n-p+1:n = λ ′ n-p+1:n , one has

where

so one shows recursively, using (8.13) and (8.16), that:

where the first term is such that

according to Lemma 11, and the second term is such that

and one concludes using (8.11) and taking

Obviously, the boundness condition on the realizations y n is not entirely satisfactory, as the generating process of (Y n ) is such that Y n should leave any interval eventually. However, Y n is marginally a Gaussian variable with variance uniformly bounded in time (since h < 1), so this remains a reasonable approximation if C y is large enough. Generalizing the above result to more general conditions is left for future research.

Application to standard state-space models

Consider a 'standard' state-space model, based on a linear auto-regressive state process (X n ): . . . i.i.d. (8.20) for t ≥ 0, ρ ∈ (-1, 1) and X 0 = Λ 0 , and an observed process (Y n ), with conditional density, with respect to an appropriate dominating measure, and conditional on X n = x n , given by the potential function Ψ X n (x n ). In this section, we show how to apply our stability results to such a standard state-space model, where the potential function depends only on the current state X n . We rewrite the model as a state space model with hidden Markov chain (Λ n ), and observed process (Y n ) corresponding to potential function

where the argument x n in the right hand side has been substituted with the appropriate function of λ 0:n , as derived from (8.20).


Guillaume Mandil
email: guillaume.mandil@usherbrooke.ca

Alain Desrochers

Alain Rivière

Framework for the Monitoring of Functional Requirements Along the Product Life Cycle

Keywords: Functional requirements, PLM, Life Cycle, Part deformation, TTRS MOTS-CLÉS : Exigences fonctionelles, PLM, Cycle de vie, Déformations, SATT

This paper presents a framework for the calculation of functional requirements at several stages of the product life cycle. It uses the TTRS/MGRE concept for the determination of dimension chains. The geometry of the assembly is represented using a parametric model constituted of points and vectors. Finally, it also integrates the calculation of part deformation due to changing operating conditions thanks to existing techniques such as finite elements.

RÉSUMÉ. Ce papier présente le cadre de développement d'une méthode pour le calcul d'une exigence fonctionnelle donnée à plusieurs étapes du cycle de vie du produit. Cette méthode utilise les concepts de SATT et EGRM pour déterminer les chaines de cotes. La géométrie est représentée grâce à un modèle paramétrique constitué de points et de vecteurs. Enfin, le calcul des déformations mécaniques dues au changement des conditions d'utilisation sera intégré par l'utilisation de techniques existantes comme les éléments finis.

Introduction

In most cases, a mechanical product is subjected to themo-mechanical loads which vary along its life cycle. These variations are inducing elastic deformations which in turn can influence the value of the functional requirements. Generally the useful value of a functional requirement is the one taken under operating conditions. However, the great majority of produts are designed and represented in Computer Aided Design (CAD) systems at the assembly stage of their life cycle. This means in fact that the values of the functional requirements are not necessarily calculated under operating conditions.

A typical application that illustrates best the above idea would be that of a jet engine for which the functional requirements varies during its own lifecycle. Indeed, the clearance between the rotor blades and engine housing (or stator) of the turbine will be quite different at assembly and in operation due to the high temperature and rotation velocity to which the rotor is subjected in service.

Therefore, a comprehensive design process would have to integrate the possibility for the designer to specify the product at the desired stage of the product life cycle. This paper proposes a framework for the calculation and the representation of a functional requirement along the product life cycle. To effectively help the designer in this way, the proposal approach integrates part deformation (as calculated from the stress analysis) into the tolerancing process.

Literature review

Analysis and synthesis of tolerance zones

In the field of tolerancing, a lot of work has been done to analyse or synthesise tolerances [Ghie et al.2007, Laperrière et al.2002] at a given phase (usually the assembly, or the design phase) of the product life cycle. Optimization of tolerances has also been carried out by several authors [Chase andGreenwood1988, Hu andPeng2007]. Some of them are optimizing tolerance in order to maximize mechanism robustness [Parkinson2000] or to minimize its cost [Singh et al.2006]. Others are considering both issues in their optimization process [Anselmetti2006].

Tolerancing deformable parts

Here are presented prior researches which have caught the authors' attention. Firstly, Samper [Samper andGiordano1998, Samper andGiordano2003] presents an approach which considers the influence of both part deformation and fit of joint into the analysis or synthesis of tolerance zones. Secondly, Cid [Cid et al.2005] developed a model which permits the evaluation of clearances under loads using a clearance torsor introduced in [Cid et al.2004]. This study investigates the case of the clearance between a vehicle door and its frame. The representation of parts considers 3D surfaces instead of 3D volumes. Finally, Pierre [Pierre et al.2008, Pierre et al.2009] has investigated part of the jet engine particular issue mentioned in the introduction.

Configuration management in PLM systems

In the field of Product Life-cycle Management (PLM), some researchers have introduced interesting concepts like Zina [Zina et al.2006] who defined the concept of "context" which could be used to define loads and environment in the proposal approach. Alternatively Eynard [Eynard et al.2004] presents an object oriented approach to help the designing team with the transmission of both design and calculation data such as geometry, use cases, loads, etc.

Product structure

In order to properly represent tolerance analysis or synthesis along the life cycle, this framework proposes the use of a specific representation for each physical state (life-cycle stage) of the product. These representations could be expressed through several views of the assembly and different representations of the parts themselves in Computer Aided Design (CAD) systems.

Assembly representations

During the design phase, different design tasks must be completed and each one uses specific tools which are in turn based on one or several representations and theories. In the perspective of the proposed framework, all the models related to geometric modeling, tolerance analysis and simulation of part deformation need to be integrated. Some of these representations are presented in figure 1.

As it is the base of product design, the nominal CAD model will occupy a central place in the proposed framework. This model translates customer's expectations (in terms of functionalities) into the geometry of a virtual product. The approach also includes the creation of an upper level of geometric representation for functional requirements such as a skeleton representing ideal functional elements such as datums, interfaces, gauges, ... These elements compose a parametric model which is used to generate the final nominal geometry. The nominal geometry (or the parametric model) is in turn exploited for other designing tasks such as tolerancing or stress analysis.

Moreover, this study integrates Finite Elements Analysis (FEA) for the calculation of part deformations due to themo-mechanical loads. Since it remains possible to obtain (or to develop) a Finite Element code for each physics involved along the product life cycle, this method can provide the designer with part deformation related to the changing environment in which the product evolves during its life-cycle. In this model, parts and assemblies are represented by discrete meshes.

Hence, the proposed framework includes several models to deal with functional requirement issues. As already explained in a previous paper [Mandil et al.2009] this framework will take functional requirements into account considering the evolution of their mean value along the product life cycle. However, it remains necessary to include tolerance zones in the model as they provide the designer with intervals for the potential value of a requirement (due to machining uncertainties). The analysis or synthesis of these intervals will be done using techniques presented in section 2. In order to monitor the evolution of the value of a given functional requirement, the framework proposes the use of TTRS 1 and MGRE 2 to build a vector loop linking the tolerance zones mentioned above. The evolution of the product along it own lifecycle will be characterized thanks to the evolution of the length and orientation of the vectors involved in the loop. The TTRS Model is based on a binary and recursive association of two functional surfaces (or group of surfaces). Globally, the goal of this association is to link each functional surface to another. The result of the association process on all the functional surfaces of a part or a mechanism is generally represented as a hierarchy tree. Additionally, this approach uses the concept of MGRE to obtain a mathematical representation of a given TTRS. The hierarchy tree is constructed by going through independent cycles in the kinematic graph. This theory is detailed in [Desrochers1991, Clement et al.1997, Desrochers2003].

Links and data exchange

In accordance with the previous subsection, and in order to properly represent assemblies, several models or views of the product will be needed. Since all of these models represent the same product, the data contained in each model has to remain consistent with the others. Therefore communication is required between models. Figure 2 shows how these models are organized as a comprehensive set of data. Among Figure 2. Links and data exchange between models the four models presented in figure 2, the parametric model supports the core of the proposed approach because the calculation of the functional requirements are made on this model. Hence, the following paragraphs will outline how this model communicates with the three others. Considering that this framework uses finite elements tools as they are, the link between the CAD model and the finite elements tools will not be detailed here.

First of all, the functional requirements will be expressed as acceptable zones, datums, interfaces, gauges, etc. which will furthermore be represented in the parametric model using points and vectors and the TTRS / MGRE approach [Desrochers1991].

Secondly, the parameters defined in the parametric model can naturally be used in the CAD software to generate the corresponding geometry. Conversly , it still remain possible to measure the geometric features of the CAD model to determine the values of some parameters.

Thirdly, after the calculation of mechanical strains through finite elements techniques, the results have to be integrated into the parametric model. To this end, the authors proposes to measure a deformed B-Rep model [Louhichi2008] and extract updated coordinates for the vectors defining the parametric model. For measurement purposes, the use of techniques inspired from metrology might prove to be helpful to determine MGRE coordinates on the deformed B-Rep model.

Product requirements across the product lifecycle

The product life-cycle defines several use cases of the product which occur under different circumstances. In a typical scenario, the product works sequentially under various conditions. As these conditions induce dimensional variations (due to variations of the thermomechanical loads applied to the mechanism), it becomes useful to obtain a representation of the product at each relevant moment of its life-cycle in order to evaluate the functional requirement under each of those conditions.

Hence, this approach proposes the use of a specific mock-up for each relevant state of the product. Those mock-ups can be viewed as different configurations of the product and managed with existing PLM 1 techniques or tools.

Moreover, since the environment is not a design variable in itself, the specifications for a given requirement under two different states have to be compatible. To ensure such a compatibility, it might be useful for the designer establish priorities between the use cases of the product. Hence, strong requirement will be expressed for high priority states wheras slack specifications will correspond to low priority ones. For the jet engine example mentioned in the introduction (section 1) a strong requirement such as a minimum and maximum gap (with close values) should be expressed at the state where the engine is running at maximum power whereas a slack requirement of noninterference might be enough to ensure assemblability. In this case, the requirement in service has a greater priority than the same one at the assembly stage.

Finally, one can notice that there exists some redundancies in the information contained in the models representing the product at several states. For example, information such as joint topology or contact surfaces are supposed to be the same under all product states. It then appears interesting to consider the use of a structure to organize the data as a compact model.

Application case

In this section, the presented framework will be illustrated on a simplified camshaft and cylinder head mechanism presented in figure 3 2 . Each step of the proposal methodology will be treated manually with the following hypothesis.

The functional requirement under scruting will be the alignment of the valve support hole (noted S3 in figure 3) with the corresponding cam. In this scenario, the parts dimensions are supposed to be known after the assembly state of the life cycle, and the designer is concerned with the value of the misalignment in service due to the thermal expansion of the two parts. The cam which is the farthest from the shoulder of the camshaft bearing has been chosen in order to tackle the most unfavorable case in terms of thermal expansion. The value of thermal expansion misalignment must remain in-1. Product Life-cycle Management 2. Enlargement in Annex I Figure 3. Cylinder head and camshaft assembly ferior to 0.2 mm under operating conditions. In our example, part deformation will not be calculated using a finite elements analysis, but rather employing a theoretical formulation.

Hypothesis

All usefull dimension values are presented in figure 3. The temperature at the assembly stage of the life cycle is 20˚C. The cylinder head is made of aluminium and the camshaft of steel. Typical values of thermal expansion coefficients are presented in table 1. In service, the cylinder head is subjected to a temperature of 90˚C whereas the camshaft temperature is 80˚C.

Material

Notation Value Aluminium α a 2.38 • 10 -5 K -1 Steel α s 1.20 • 10 -5 K -1
Table 1. Typical values of thermal expansion coefficient

In order to keep calculations simple, the authors have considered the restrictive but reasonable hypothesis of a one dimensional thermal expansion of the two parts along the cam shaft axis. The dimensions after thermal expansion will be computed thanks to the equation 1 where S f and S i stand for final and initial stage of the assembly life cycle, T for temperature and L for length. One can note that temperature and length depend on the life cycle stage.

L(S f ) -L(S i ) = α • L(S i ) • (T (S f ) -T (S i ))
[1]

Parameterization of the model

The first step in obtaining the parameterization vectors is to associate a TTRS to each pair of functional surfaces. Secondly, one has to build the TTRS tree for each part using the association rules described in [Desrochers1991, Clement et al.1997, Desrochers2003]. The construction of the TTRS tree induces the choice of an associated MGRE for each TTRS. Finally the vectors are deduced from both the TTRS tree and their associated MGRE.

For the construction of the TTRS tree representing the camshaft as shown in figure 4a , one must first denote that all the surfaces involved are planar. Consequently surfaces are associated as planar TTRS and their corresponding MGRE are matching theoretical planes. More specifically the example features three planar surfaces : B, S1 and S2. In the TTRS tree construction process, S1 and S2 have been associated first. As these two planes are parallel, the TTRS resulting from this association (P4) is also a planar TTRS parallel to its children surfaces. Accordingly, the associated MGRE is the median plane of S1 and S2 bears the name P4. Secondly, surface B is associated with P4, and the resulting TTRS (called CS in figure 4) is also planar and parallel to B and P4. The TTRS CS carries the distance between shoulder B and MGRE P4. Finally, the associated vector in the parametric model is normal to B and its norm is equal to the distance between B and P4 ; it is noted --→ V CS .

The same method is used to build the TTRS tree of the Cylinder head (figure 4b). This part has two functional surfaces, plane A and Cylinder S3. Cylinder S3 is represented by its axis or MGRE as the A5 line. Plane A and A5 are associated as TTRS CH. As the A5 line is parallel to plane A, the resulting TTRS is prismatic. This TTRS obviously carries the distance between plane A and A5. Therefore, the associated vector in the parametric model is normal to plane A and its norm is equal to the distance between plane A and line A5 ; it is noted --→ V CH .

This parameterization allows the designer to characterize the alignment default between the cam and the valve support hole using to the norm :

--→ V CS ---→ V CH

Calculations

Thanks to the hypothesis in section 5.1 and the dimensions specified in figure 3, it becomes possible to calculate with equation 1, the norm of --→ V CS and --→ V CH under their service temperatures. The results are presented in table 2. They show that during from the assembly stage of the life cycle to operating conditions the cam and the valve support hole are subjected to unacceptable misalignment. Indeed this default is greater than 0.2 mm.

As mentionned in section 4, the designer must establish priorities between the several stages of the life cycle. In this case, it appear obvious that the respect of functional requirement under operating conditions is more important than at assembly. Therefore it becomes interesting to compute what dimensions must be specified at 20˚C in or- These results (in table 3) show that there is no misalignment anymore under operating condition. However the 0.235 mm wide misalignment default is now detectable at assembly state (at 20˚C). This misalignment default has to be added to the default calculated during the analysis of the functional requirement along the corresponding dimension chain (which has not been carried out in this study). If the accumulation of these two sources of misalignment produces a non acceptable default at assembly stage, then the complete system has to be redesigned altogether. In order to reduce

  1. Technologically and Topologically Related Surfaces 2. Minimum Geometric Reference Element 3

Figure 1 .

 1 Figure 1. Models used for the representation of the mechanism

Figure 4 .

 4 Figure 4. TTRS tree of the parts Vector Final temperature Norm at 20˚C Norm at final temperature --→ V CS 80˚C 375 mm 375.270 mm --→ V CH 90˚C 375 mm 375.625 mm --→ V CS ---→ V CH 0 mm 0.355 mm

Table 2 .

 2 Final norm of parameterization vectorsder to have the cam and the valve support hole aligned under operating conditions. To compute these dimensions, the same approach was used as that described earlier. Hypothesis and results are presented in table 3.

	Vector Initial temperature Norm at initial temperature Norm at 20˚C
	--→ V CS	80˚C	375.625 mm	375,355 mm
	--→ V CH	90˚C	375.625 mm	375 mm
	--→ V CS -	--→ V CH	0 mm	0.355 mm

Table 3 .

 3 Updated norm of parameterization vectors at 20˚C

 CPI'2009.the misalignment due to thermal expansion, one solution would be to place the shaft shoulder in the middle of the camshaft itself, therefore distributing the thermal expansion on both sides.

Conclusion and perspectives

This paper has presented a framework for the representation and calculation of functional requirement values along the product life cycle. The proposed framework is based on existing concepts used in the field of tolerancing such as TTRS and MGRE for the definition of dimension chains and the study of functional requirements. It also uses existing tools like finite elements for the calculation of part deformation. It finally proposes the original idea of representing, with distinct configuration, each of the product during its life cycle. To this end, a point and vector based model originating from the TTRS/MGRE representation is build. The advantage of this approach stems from the use of high level management tools such as TTRS and MGRE for geometries and tolerance zones representation. The future development of this approach will include the integration of 3D dimension chains and the automation of the method. The authors also envision the integration of this framework with PLM tools for the management of configurations, and the use of an object oriented formalism to represent parts and mechanisms and organize the data contained at the different levels of product modelisation.

Annex I : Enlargement of figure 3