J Juillard
email: jerome.juillard@supelec.fr

A Bonnoit

E Avignon

S Hentz

E Colinet

Large amplitude dynamics of micro/nanomechanical resonators actuated with electrostatic pulses

In the field of resonant NEMS design, it is a common misconception that large-amplitude motion, and thus large signal-to-noise ratio, can only be achieved at the risk of oscillator instability. In the present paper, we show that very simple closed-loop control schemes can be used to achieve stable largeamplitude motion of a resonant structure, even when jump resonance (caused by electrostatic softening or Duffing hardening) is present in its frequency response. We focus on the case of a resonant accelerometer sensing cell, consisting in a nonlinear clamped-clamped beam with electrostatic actuation and detection, maintained in an oscillation state with pulses of electrostatic force that are delivered whenever the detected signal (the position of the beam) crosses zero. We show that the proposed feedback scheme ensures the stability of the motion of the beam much beyond the critical Duffing amplitude and that, if the parameters of the beam are correctly chosen, one can achieve almost full-gap travel range without incurring electrostatic pull-in. These results are illustrated and validated with transient simulations of the nonlinear closed-loop system.

Introduction

Resonant sensing consists in measuring the frequency shift of a system subject to the variation of a given physical quantity. Because of its moderate complexity, this measurement technique is becoming commonplace in the context of MEMS and NEMS devices [START_REF] Stemme | Resonant silicon sensors[END_REF][START_REF] Roessig | Integrated MEMS tuning fork oscillators for sensor applications[END_REF][START_REF] Ekinci | Ultimate limits to inertial mass sensing based upon nanoelectromechanical systems[END_REF]. This paper focuses on closed-loop resonant sensors, where the micromechanical structure is brought to oscillate by being placed inside a feedback loop [START_REF] Gallacher | A control scheme for a MEMS electrostatic resonant gyroscope excited using combined parametric excitation and harmonic forcing[END_REF][START_REF] Sung | Design and performance test of an oscillation loop for a MEMS resonant accelerometer[END_REF][START_REF] Ahn | Capacitive microbeam resonator design[END_REF]. In the present work, the structure is a clamped-clamped beam, the motion of which is sensed capacitively. It is maintained in an oscillation state with pulses of electrostatic force that are delivered whenever the detected signal (the position of the beam) crosses zero [START_REF] Colinet | Actuation of resonant MEMS using short pulsed forces[END_REF][START_REF] Juillard | From MEMS to NEMS: closed-loop actuation of resonant beams beyond the critical Duffing amplitude[END_REF]. In order to maximize the signal-to-noise ratio (SNR) and, thus, to relax the constraints on the electronic design, the detected signal must be as large as possible, which means that, for a given set of structural parameters and a given bias voltage, the oscillation amplitude of the resonant beam must also be as large as possible. This raises questions concerning what oscillation amplitude can be sustained without incurring mechanical [START_REF] Kaajakari | Nonlinear limits for single-crystal silicon microresonators[END_REF][START_REF] Lifshitz | Nonlinear dynamics of nanomechanical and micromechanical resonators[END_REF] or electrostatic [START_REF] Leus | On the Dynamic Response of Electrostatic MEMS Switches[END_REF][START_REF] Nielson | Dynamic pull-in of parallel-plate and torsional electrostatic MEMS actuators[END_REF] instability. We show in this paper that, in spite of the nonlinearities, the proposed feedback scheme ensures the stability of the motion of the beam much beyond the critical Duffing amplitude. We also show that, if the parameters of the beam are correctly chosen, one can realistically achieve almost full-gap travel range without incurring electrostatic pull-in. In section 2, the sensing cell of the resonant accelerometer that is the basis of our work is briefly described. The focus is brought on the capacitive detection and actuation schemes. In section 3, a simplified one-degree-offreedom (1-DOF) model of the beam is derived. In particular, an approximate expression of the first modal component of the electrostatic force acting on the beam is determined. In section 4, the electrostatic instabilities that are inherent to the design are addressed. We show, using the model of section 3, that it is possible to achieve almost full-gap travel range without pull-in if the bias voltage is correctly chosen. In section 5, describing function analysis (DFA) is used to characterize the oscillations of the closed-loop system. We show that, when the electrostatic pulses are short with respect to the natural period of the system, there exists a unique oscillation regime and that it is stable, provided the condition for electrostatic stability established in section 3 is met. These results are illustrated in section 6 and validated with transient simulations of the closed-loop system.

Framework

Sensing principle

The accelerometer consists in a suspended seismic mass, which is anchored to the substrate via a short beam acting as a flexure. A longer, more flexible beam, perpendicular to the flexure, acts as a sensing element (figure 1). When the mass is subject to acceleration, the sensing beam undergoes compressive or tensile stress (depending on the direction of the acceleration), resulting in a variation of its stiffness. This variation can then be detected via a shift ω ∆ of its natural pulsation 0 ω .

Two sensing principles may be used to monitor frequency shifts. In open-loop sensing, the structure is excited using a repeated chirp (a sine wave with slowly increasing or decreasing frequency). The frequency response of the structure can then be calculated and the value of the "open-loop resonance frequency" (the frequency for which the frequency response is maximal) can then be determined and monitored. In closed-loop sensing, the micromechanical structure is placed inside a feedback loop, without external input. If the feedback loop is correctly designed, the structure starts oscillating at a so-called "closed-loop resonance frequency", which shifts along with the natural frequency of the structure.

Description of the accelerometer sensing cell

In the ANR-funded M&NEMS project, which sets the framework of this paper [START_REF] Juillard | From MEMS to NEMS: closed-loop actuation of resonant beams beyond the critical Duffing amplitude[END_REF], the sensing beam is placed between two equidistant electrodes, so as to be the midpoint of a capacitive half-bridge. The transmission line and contact pads between the structure and the control electronics give rise to a parasitic capacitance p C , typically on the order of several pF, which is much larger than the nominal capacitance of the structure (on the order of 1fF) and leads to poor SNR [START_REF] Colinet | Ultra-Sensitive Capacitive Detection Based on SGMOSFET Compatible With Front-End CMOS Process[END_REF]. At the output of the charge amplifier, the voltage out V can be written:

offset e f p out V Q C K C K V + + + = ) 1 ( , (1)
where e Q is the charge accumulated on the mobile beam, K is the gain of the operational amplifier, f C is the feedback capacitance and offset V is an integration constant that is removed with high-pass filtering. We may then assume throughout the rest of the paper that 0

= offset V and that 0 = e Q
when the beam is halfway between the electrodes, so that 0 = out V when the displacement of the beam w is uniformly 0. In order to maximize the SNR, the detected signal must be as large as possible, which means that, for a given set of structural parameters and a given bias voltage, the oscillation amplitude of the resonant beam must also be as large as possible. We show in what follows that such large-amplitude motion may be sustained without incurring mechanical or electrostatic instability, if a proper feedback scheme is used to actuate the beam.

Actuation and detection schemes

Applying brief negative or positive pulses of electrostatic force to the beam whenever it passes through a certain reference position (typically

0 = = w V out
) can bring the system into a sustained self-oscillation state [START_REF] Colinet | Actuation of resonant MEMS using short pulsed forces[END_REF][START_REF] Juillard | From MEMS to NEMS: closed-loop actuation of resonant beams beyond the critical Duffing amplitude[END_REF]. The voltage pulses (with amplitude p V ± and duration

0 2 ω π << p T
) can be delivered to the central electrode, as in [START_REF] Colinet | Actuation of resonant MEMS using short pulsed forces[END_REF], or through the biasing electrodes (figure 1). Both solutions are equivalent from a dynamical system point of view. In either case, provided the length L and width b of the beam are large with respect to G , the electrostatic gap, one may consider that every slice of area bdx of the beam is subject to an electrostatic force dx F e , where lineic density e F is given by the plane capacitance approximation:

( ) ( ) ( ) ( ) ( ) ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ - ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ - - = ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + + - - - = w V G V G V w V w G b w G V V w G V V b F c b c b c b c b e 2 2 2 2 / 2 / 2 2 2 2 0 2 2 2 2 0 ε ε , (2)
where 0 ε is the permittivity of vacuum, b

V is the bias voltage and c V is the control voltage (i.e. the voltage applied to the mobile beam). Under the same assumptions, the lineic charge density accumulated on an elementary slice of the beam is:

⎟ ⎠ ⎞ ⎜ ⎝ ⎛ - - - = ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + + + - + - = 2 2 2 2 0 0 2 2 / 2 / w G w V w G G V b w G V V w G V V b b c c b c b e ε ε σ , (3)
When no pulse is being delivered to the beam, 0 = c V and (2) and (3) simplify to:

( ) 2 2 2 2 0 2 w G w bGV F F b soft e - = = ε , (4)
which corresponds to the electrostatic softening of the system, and (

) 5
Note that (4) and (5) vanish when w is uniformly 0 or, assuming that the deformation of the beam is unimodal, when 0 = out V .

Since the pulses are very short, there is only a very small amount of crosstalk between the actuation and detection schemes. Assimilating these finite pulses to ideal Dirac pulses, one may then consider that (4) and (5) are always valid except at instants 0 t when out V is exactly zero, in which case: . [START_REF] Colinet | Actuation of resonant MEMS using short pulsed forces[END_REF] The fact that, in practice, the pulses have finite duration is not detrimental to the proper functioning of the system, because when the self-oscillation conditions are met, voltage pulses are delivered in such a way that they do not modify the sign of out V [START_REF] Colinet | Actuation of resonant MEMS using short pulsed forces[END_REF]. In other words, the sign of out V is always the exact opposite as the sign of the displacement of the beam, even during the voltage pulses. This makes the design of the electronic feedback loop very simple: basically, all it takes to generate voltage pulses when 0 = w is a comparator and a pair of monostable multivibrators.

p p b t t e T V V G b dt F 2 0 0 0 0 lim ε τ τ τ ± = ∫ + - → (6)
Assuming the electrodes have the same length and width as the beam, the system is governed by: ( )

pulse b F t w bh t w w G w bGV x w T x w Ebh = ∂ ∂ + ∂ ∂ + - - ∂ ∂ - ∂ ∂ 2 2 2 2 2 2 0 2 2 4 4 3 2 12 ρ µ ε (8),
where E is Young's modulus, ρ is the density, h is the thickness of the beam, µ is a linear damping coefficient, T is the axial force caused by the elongation of the beam:

∫ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ∂ ∂ = L dx x w L Ebh T 0 2 2 . ( 9
) soft e pulse F F F - =
is the "pulsed" component of the electrostatic force and consists of Dirac pulses occurring

when 0 = out V
, with amplitudes given by ( 6). In the next section, we define the pull-in voltage of the structure and its pull-in amplitude and show that, with proper design, one can use pulse-actuation to achieve almost full-gap travel range.

Model derivation

In the context of a resonant sensor application, it is a reasonable to assume that the deformation of the beam is unimodal. Let us then write ( ) ( ) ( )

x w t a t x w , ~0 =
, where 0 w is the first clamped-clamped beam eigenmode, with

eigenvalue 4 0 β ( 730 . 4 0 ≈ β
) and where the following dimensionless quantities are used:

G w w / ~= , L x x / ~= and t t 0 ~ω = , ρ ρ β ω E L h E L h 2 2 2 0 0 459 . 6 12 ≈ = . ( 10
)
The value of ( )

2 / 1 0 w is set to 1, so that 1 ± = a
corresponds to full-gap displacement of the midpoint of the beam. Projecting (8) on 0 w , we obtain the following model of the beam dynamics:

( ) [ ] ( ) a a f a Q a a I a a , 1 2 ∆ = + + - + δ γ , (11)
1 The terms involving V in (2) are cancelled out because the pulses are delivered when 0 = w . Thus, the amplitude of the

actuation is proportional to p V , regardless of whether b p V V << .
where ( )

a a & , ∆
designates a Dirac pulse occurring when 0 = a , with the same sign as a & (giving the oscillator a kick in the direction of its motion),

2 2 1 0 2 0 2 1 0 2 0 4 0 2 2 719 . 0 6 h G x d w x d x w h G ≈ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ ∂ ∂ = ∫ ∫ β γ , (12) 2 2 pi b V V = δ , 2 / 1 4 0 3 3 2 / 1 4 0 3 3 2 0 134 . 9 6 ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ≈ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ = L G Eh L G Eh V pi ε ε β , (13) ( ) ( ) ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - = ∫ ∫ 1 0 2 0 1 0 2 2 0 2 2 0 1 x d w x d w a w a I , (14)
p pi p b p pi p b p V V V x d w x d w V V V f τ τ 2 1 0 2 0 1 0 0 2 639 . 2 2 ≈ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ = ∫ ∫ , (15)
where

p p T 0 ω τ =
and Q is the quality factor of the system. There exists no closed-form expression of ( ) a I , however it is possible to approximate it with:

( ) ( ) ( ) 2 / 3 2 1 ˆ-- = ≈ a a I a I , (16)
as is explained in Appendix A. The relative error between ( ) a I and ( )

a I ˆ
is less than 5% across the whole gap -it can even be reduced to less than 1% by multiplying ( )

a I ˆ with ( ) 2 041 . 0 1 a - (see Appendix A).
Because of its range of validity, this approach is better suited to large amplitude dynamics than Taylor series expansion, which is commonly used in MEMS/NEMS modeling. In particular, it allows us to capture the dynamic pull-in behavior of a pulse-actuated in a simple way, as is shown in section 4.

Pull-in phenomena

A fundamental limit to the displacement amplitude of electrostatically-actuated devices is set by the static pull-in phenomenon [START_REF] Leus | On the Dynamic Response of Electrostatic MEMS Switches[END_REF][START_REF] Nielson | Dynamic pull-in of parallel-plate and torsional electrostatic MEMS actuators[END_REF]. The pull-in amplitude is classically determined as the displacement amplitude for which the electrostatic softening force and the mechanical restoring force are balanced but the equilibrium is unstable. In the case of the clamped-clamped beam setup of figure 1, it is clear that, in the absence of (pulsed) actuation forces, 0 = w is always an equilibrium point and that, provided b V is small enough, this equilibrium point is stable. The threshold value of b V beyond which this equilibrium becomes unstable is determined in sub-section 4.1. We then show in sub-section 4.2 that the electrostatic nonlinearity sets an upper bound to the oscillation amplitude of the structure and determine this limit.

Determination of the static pull-in voltage

The problem of the beam at rest subject to an increasing bias voltage is similar to the problem of a beam subject to increasing compressive stress: the beam is gradually softened until a critical point is reached for which it is pulled-in (or it buckles, in the mechanical analogy). As in classic Euler buckling problems [START_REF] Chow | Methods of Bifurcation Theory[END_REF], the equilibrium position ( 0 = w ) becomes unstable when the first eigenvalue of the linearized spatial partial differential operator appearing in (8) equals zero. In other words, the beam is pulled-in from its central position when the linearized modal stiffness of the system equals zero. From [START_REF] Leus | On the Dynamic Response of Electrostatic MEMS Switches[END_REF], the modal stiffness is ( )

a I a δ γ - + 2 1
. Linearizing this expression for 0 = a yields the following pull-in condition:

1 = δ , (17) i.e. pi b V V ± = .
Working with 1 < δ ensures only that 0 = w is a stable equilibrium point around which the beam can oscillate. In the following sub-section, we show that there exists a maximal displacement amplitude beyond which the beam is pulled-in to one of the electrodes and that this amplitude depends on b V .

Determination of the maximal displacement amplitude

The first term on the left-hand side of (11) corresponds to the restoring forces acting on the beam, including the linear elastic and nonlinear elastic and electrostatic components, which derive from a potential energy pot E .

Using [START_REF] Kacem | Nonlinear dynamics of nanomechanical beam resonators: improving the performance of NEMS-based sensors[END_REF], pot E can be expressed as:

2 4 2 1 1 4 1 2 1 a a a E pot - - + = δ γ . ( 18
)
pot E is represented in figure 2 for different values of δ and γ .

When 1 < δ

, the central position is stable, as shown in sub-section 4.1. pot E is maximal when:

( ) 0 1 1 2 / 3 2 2 = - - + - a a δ γ , ( 19
) which has two solutions max a ± ( 0 max > a
), both of which correspond to unstable equilibria. Now, let us suppose that energy is injected into the system (initially at rest) whenever 0 = a , in the form of short electrostatic pulses, so that an oscillation starts building up. If too much energy is injected into the system, it moves past the potential barrier in max a ± and the beam is inevitably pulled-in to one of the electrodes. Otherwise, it eventually reaches a steady-state regime and stabilizes at an amplitude which is inferior to max a . Thus max a is the maximal displacement amplitude that can be achieved by an oscillating pulse-actuated clamped-clamped beam without incurring pull-in. When 0 = γ , the following analytical expression can be derived from (19):

3 / 2 max 1 δ - = a . ( 20
)
Accounting for stiffening results in a larger value of max a than predicted by (20) (figure 3), however no analytical expression of max a can be found. It is notable that very small values of δ result in max a being very close to 1. 1 ≥ δ corresponds to the case when the central position is unstable (figure 2). However, the fact that the beam is "pulled-in" does not mean it necessarily comes in contact with an electrode: it might instead get stuck in one of the potential wells that appear for large values of γ . In this configuration, it should also possible to use a pulse- actuation scheme to achieve stable oscillations around 0 = a even though the central position is intrinsically unstable.

In the following section, we assume that 1 < δ .

Analysis of the closed-loop system

Describing function analysis, or equivalent linearization, is a method used (mostly) by engineers to determine certain properties of weakly nonlinear systems subject to various stimuli. When the only stimuli are sine waves, the method is also known as harmonic balance, in which case it yields the same results as those one would obtain through the method of averaging (or KBM method) [START_REF] Gelb | Multiple-input describing functions and nonlinear system design[END_REF] used to first order. In this context, the difference between the two methods is purely a matter of perspective, the describing function formalism being more oriented toward the design of control systems. 5.1. Describing function analysis of the closed-loop system An equivalent block-diagram representation of (11) is shown in figure 4. One may consider the closed-loop system as a linear system with a triple feedback nonlinearity. If the "filter hypothesis" holds (i.e. if the linear part efficiently filters out the high-order harmonic content output by the nonlinearity), one may apply DFA to the system and expect good qualitative and quantitative results [START_REF] Gelb | Multiple-input describing functions and nonlinear system design[END_REF]. DFA is then particularly well-suited to the case of a resonant structure with a large quality factor oscillating with a pulsation close to 0 ω .

The existence of a periodic regime, and its characteristics (amplitude, pulsation and stability) can then be determined by: 1. assuming

( ) t A a ω sin = , 1 0 < < A , 0 > ω 2. finding the equivalent complex gain ( ) ω , A N of the nonlinear part ( ) a a F & , , ( ) ( ) ( )
∫ = ω π ω π ω ω / 2 0 exp , , dt t j a a F A A N & (21)
3. checking for an oscillatory regime satisfying the Barkhausen condition, i.e. solving for the limit cycle amplitude osc A and pulsation osc ω ( ) ( )

0 1 , = + ω ω j H A N (22)
where ( )

ω j H
is the transfer function of the linear part. (22) can be split in two equations and rewritten as

( ) ( ) 0 , Im , Re = = ω ω A A ( 23
)
where Re and Im are the real and imaginary part of ( ) ( )

ω ω j H A N / 1 , + .
Step 3 is equivalent to imposing that the total phase lag in the loop is equal to zero. An oscillatory regime defined by (22-23) is stable with respect to slight perturbations if the following inequality is satisfied:

0 Re Im Im Re > ∂ ∂ ∂ ∂ - ∂ ∂ ∂ ∂ ω ω A A , ( 24
)
where the partial derivatives are taken at The nonlinear feedback block consists in three nonlinearities in parallel. Their equivalent complex gain is then given by: ( )

( ) ( ) ( ) ( ) ω ω , , A N A N A N A N pulse soft Duff + + = , (25) where ( )
2 4 3 A A N Duff γ = , (26) and
( )

p pulse f A j A N π ω ω 2 , - = . ( 27
)
One can directly compute an approximation of ( )

A N soft
, without passing through the intermediate step of (16) 2 .

This yields:

( ) ( ) ( ) ( ) ( ) 2 2 / 2 0 1 1 exp sin sin , A A dt t j t t A I A N soft - + - ≈ - = ∫ κ δ ω ω ω π ω δ ω ω π , (28)
where 182 . 0 ≈ κ (cf. appendix A). The Barkhausen criterion (22) then becomes:

Q j f A j A A A p ω ω π ω κ δ γ - - = - - + - 1 2 1 1 4 3 2 2 2 2 . ( 29
)
There exists one unique couple ( )

osc osc A ω ,
that cancels out the real part and the imaginary part of (29):

2 There is little point in using ( )

a I ˆ
to derive an analytical expression of ( ) A N soft : this approach leads to an expression involving elliptic integrals which can only be apprehended with difficulty (and should then be approximated, leading to a simpler but less accurate expression which would be an approximation of an approximation).

⎪ ⎪ ⎩ ⎪ ⎪ ⎨ ⎧ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - + - + = = 2 / 1 2 2 2 1 1 4 3 1 2 osc osc osc osc p osc A A A Q f A κ δ γ ω π . (30)
The values that may be taken by osc A are theoretically limited by the second equation, i.e. there can be no oscillation regime such that max A A osc >

, where max A satisfies:

0 1 1 4 3 1 2 max 2 max 2 max = - + - + A A A κ δ γ , (31)
for which 0 = osc ω . However, one should not expect this limit to have a practical importance: first of all, in many practical cases max max A a < : this is typically the case when 0 > γ (hardening Duffing behavior) and 1 0 < < δ . Thus, the electrostatic instability discussed in section 3.2 sets a tighter limit on the oscillation amplitude than (31). Moreover, as we have already mentioned, DFA accurately predicts the behavior of systems that satisfy the filter hypothesis: it is clear that, if osc ω decreases sufficiently, the high-order harmonics generated by the different nonlinearities are no longer filtered out by the linear part and the filter hypothesis holds no more. Thus, the value of max A that can be inferred from (31) is out of the range of validity of DFA. This shows that one should be very careful when trying to determine the onset of electrostatic instability under sinusoidal motion assumptions [START_REF] Kacem | Nonlinear dynamics of nanomechanical beam resonators: improving the performance of NEMS-based sensors[END_REF][START_REF] Fargas-Marques | Resonant pull-in condition in parallel-plate electrostatic actuators[END_REF], because the two notions (instability and sinusoidal motion) are somehow contradictory. This is confirmed by a stability analysis based on DFA: the oscillatory regime defined by ( 30) is stable provided (24) is satisfied. We have:

0 Im = ∂ ∂ osc osc A ω ω , osc A osc osc ω ω ω 2 Re - = ∂ ∂ and Q A A osc osc A osc osc ω ω = ∂ ∂ Im , (32)
thus the stability condition boils down to:

0 0 2 2 > ⇔ > osc osc osc Q A ω ω , (33)
which is always satisfied. However, we have shown in section 3 that the amplitude of the oscillations of a pulseactuated structure cannot exceed max a . Thus the resonant stability condition is simply:

max a A osc ≤ . ( 34
) 35
This value is "optimal" in the sense that small variations of the amplitude will not affect the value of osc ω : thus, one may achieve large amplitude oscillations and good frequency stability, both of which would be very desirable features in MEMS/NEMS reference oscillator applications [START_REF] Juillard | Enhancing the frequency stability of a NEMS oscillator with electrostatic and mechanical nonlinearities[END_REF][START_REF] Kaajakari | Phase noise in capacitively coupled micromechanical oscillators[END_REF].

Describing function analysis of the open-loop system

Before moving on to the comparison of these results with simulated data, we stress the fact that the existence and the uniqueness of a stable periodic regime for an electrostatic pulse-actuated clamped-clamped beam are only valid when the system is operated in closed-loop, and in particular when the pulses are triggered when the detected signal crosses zero.

A simplified model of the open-loop beam is represented in figure 5 3 . In this configuration, there may exist up to five oscillation states at a given excitation amplitude and pulsation: this can be inferred from the DFA of [START_REF] Leus | On the Dynamic Response of Electrostatic MEMS Switches[END_REF], where one assumes that ( )

ϕ ω + = t A a sin
(the fact that the loop is open no longer ensures that a is in phase with the excitation) and the origin of time coincides with a positive pulse. Projecting this equation on ( ) t jω exp yields, after using (28):

( ) ϕ ϕ π ω ω ω κ δ γ cos sin 2 1 1 4 3 1 2 2 2 2 j f A Q j A A A + = ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ + ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - - + - + , ( 36
)
where f is the amplitude of the Dirac pulses injected into the system.

Taking the squared modulus of (36) results, after multiplication by ( )

2 2 1 A - , in a fifth-order polynomial in 2 A : ( ) ( ) ( ) ( ) ( ) ( ) 0 4 1 8 2 2 1 2 4 3 1 2 1 2 2 4 4 3 1 1 2 4 3 1 2 1 2 3 1 1 2 2 3 16 9 4 3 1 2 3 16 9 2 2 2 2 4 2 2 2 2 2 4 4 2 2 2 2 6 4 2 2 2 2 8 2 10 2 = - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + + - + - + ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ + - - - - + - ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ + - - + ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ + ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ + - - + + + + - - - + ⎟ ⎠ ⎞ ⎜ ⎝ ⎛ - + - + π ω ω ω π δ δ ω ω δκ γ δ π δκ γ δ ω ω δκ γ γ δκ κ δ γ γ ω δκ γ γ γ f A Q f A Q f A Q A A . ( 37
)
This polynomial may have multiple roots, meaning that the open-loop system may have up to five oscillation states 4 and exhibit hysteretic behavior. The corresponding values of ϕ are given by:

⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - + - + = 2 2 2 2 1 1 4 3 1 tan ω κ δ γ ω ϕ A A A Q . ( 38
)
Closing the loop simply ensures that ϕ takes a specific value (in our case, 0 = ϕ ) and, thus, results in a single possible oscillation state (figure 6). The case when the sensor is operated in closed-loop in the presence of phase-delay (consequent to the filtering of the detected signal) between the pulses and the triggering event is the subject of forthcoming work. The fundamental difference between these cases and the zero-phase delay case covered in the present work is that the amplitude of the pulses is modulated by the electrostatic nonlinearity: for large values of the phase-delay, this may result in multiple oscillation states, even though the loop is closed.

Simulation and results

The purpose of this section is to validate with simulations the theoretical results obtained in the previous sections. This validation is conducted in two steps: first, the validity of the describing function approach is demonstrated by simulating the one-degree-of-freedom model [START_REF] Leus | On the Dynamic Response of Electrostatic MEMS Switches[END_REF], in which ( )

a I is replaced by ( ) a I ˆ (16
) and the pulses have finite duration. A more complex model, obtained by projecting (8) on three eigenmodes, is then simulated to verify the well-foundedness of the simpler model [START_REF] Leus | On the Dynamic Response of Electrostatic MEMS Switches[END_REF]. The parameters of the resonator correspond to those of the accelerometer structure (figure 7) developed in the ANR-funded M&NEMS project:

L =25µm, b =500nm, G h = =250nm, Q =6000, ρ =2320kg.m -3 , E =149GPa. For this set of parameters, 0 ω ≈2.07×10 7 rad.s -1 ,
3 It is simplified in the sense that, in the system represented in figure 1, the amplitude of the pulses of electrostatic force depends on the value of the phase-delay, as can be inferred from (2). 4 When 0 = δ , up to three oscillation states may be observed. This corresponds to the well-known 'critical amplitude' phenomenon associated with the Duffing pendulum [START_REF] Kaajakari | Nonlinear limits for single-crystal silicon microresonators[END_REF]. pi V ≈29.63V, γ ≈0.719. In the actual circuit, the duration of the impulses is set to p T =10ns ( 207 . 0 = p τ ), so that they may be considered short with respect to the period of oscillation. The bias voltage that is applied to the structure is equal to the maximum supply voltage of the technology chosen for the circuit -in the present case, ST CMOS 130 nm-i.e. b V =1.2V, which corresponds to δ ≈1.64×10 -3 .

Solving [START_REF] Kaajakari | Phase noise in capacitively coupled micromechanical oscillators[END_REF], we find that the pull-in amplitude is 995 . 0 max ≈ a . From (30) and (15), the beam should not be pulled-in, provided when the steady-state regime has supposedly been reached. In theory, the amplitude and the pulsation should behave according to (30), as depicted in figure 8. The relative error between simulation and theory is represented in figure 9. As anticipated, the error increases with p τ . However, the agreement is very good, even for "long" pulses, with relative errors inferior to 4% for the amplitude and inferior to 2% for the pulsation, even when the oscillation amplitude is very close to 1. This validates the describing function approach of section 5. Moreover, we verify that the beam pulls-in when the oscillation amplitude becomes larger than 995 . 0 max ≈ a .

Simulation of the 3-DOF model

A less idealized model is simulated: (8) is decomposed and projected on the first three relevant eigenmodes of the beam. The projection integrals are evaluated at each integration step with a first-order Newton-Cotes quadrature scheme. p τ -long pulses are triggered when the total charge accumulated on the beam crosses zero. The charge is calculated at every time-step by integrating (5) along the length of the beam. The oscillation amplitude is defined as the maximal displacement amplitude of the midpoint of the beam and the oscillation period is calculated as the difference between two consecutive zero-crossings of the total charge. The errors have the same orders of magnitude as for the 1-DOF model and the same behaviour with respect to p τ . This validates the use of (30) to predict the behaviour of the actual system. However, some qualitative differences do exist between the 1-DOF case and the 3-DOF case, as illustrated in figure 10. Furthermore, the beam is pulled in at 982 . . Thus, for large oscillation amplitudes, the second and third modes have a small (but finite) contribution to the elastic energy of the system: it is difficult to derive general results from these simulations without resorting to finer, more complex, models than those developed in section 5. For example, using DFA makes it impossible to predict sub-harmonic resonances that are likely to appear in lowfrequency open-loop excitation contexts (more appropriate methods are described in [START_REF] Nayfeh | Secondary resonances of electrically actuated resonant microsensors[END_REF]). The case of inter-mode crosstalk may be addressed if simplifying assumptions are made (for example, if all the modes oscillate in phase).

The system was simulated for several other sets of parameters, in closed-loop as well as in open-loop, resulting in very small errors between theory and simulation, as long as the hypotheses made in section 5 were respected.

Conclusion

In this paper, we have shown that very simple closed-loop control schemes can be used to achieve stable largeamplitude motion of a resonant structure, even when jump resonance (caused by electrostatic softening or Duffing hardening) is present in its open-loop frequency response. The case of a capacitive pulse-actuated resonant accelerometer sensing cell was thoroughly investigated. The nonlinear softening term was approximated with simple expressions, valid across the whole gap. This helped us establish the maximal displacement amplitude max a that can be achieved by a closed-loop pulse-actuated clamped-clamped beam without incurring electrostatic instability: it was shown that, depending on the structure's design parameters γ and δ , almost full-gap travel range is possible. The same approach may also be used for other beam-electrode configurations or for other actuation schemes. Using describing function analysis, we studied closed-loop and open-loop pulse-actuation and showed that the closed-loop approach ensures the stability of the motion of the beam even when the open-loop frequency response has a hysteretic characteristic. Some analytical expressions of the oscillation amplitude and of the oscillation pulsation were established and validated with transient simulations of the nonlinear closed-loop system, in the case of the sensing cell of the ANR-funded M&NEMS project.

The most immediate practical consequence of the present work is that, with proper mechanical and feedback design, very large oscillation amplitudes and, thus, signal-to-noise ratios may be achieved, thereby relaxing the constraints on the design of the analog front-end of the electronic circuitry. Another practical consequence is that, as mentioned in section 5, the softening nonlinearity can be used to compensate for the hardening nonlinearity to achieve large oscillation amplitudes together with good frequency stability [START_REF] Juillard | Enhancing the frequency stability of a NEMS oscillator with electrostatic and mechanical nonlinearities[END_REF]. The analysis of the closed-loop system when phase-delay is present in the feedback loop is the subject of ongoing work. Furthermore, an ASIC implementing the electrostatic actuation and detection scheme described in this paper has been developed and should soon be tested.

Appendix A Equation ( 16) is obtained by fitting the projection of the electrostatic force on the first eigenmode with a two parameter model parameterized by 0 n and 0 κ . Let us then define

{ } ( ) ( ) ( ) ∫ - ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - - + = ε κ ε ε κ κ 1 0 2 2 2 , 1 1 min arg , ˆda a I a a n n n , (A-1)
where ( ) ( ) The same approach can be used to determine an approximation of ( ) , which minimizes (A-1), one may find a closed-form expression (involving complete elliptic integrals, as mentioned in section 4) of the corresponding describing function gain. In the neighborhood of 1 = a , this gain behaves as: in the same region, as is indeed the case (figure A2).

⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ ⎟ ⎟ ⎠ ⎞ ⎜ ⎜ ⎝ ⎛ - = ∫ ∫ 1 0 2 0 1 0 2 2 0 2 2 0 1 x d w x d w a w a I , (
A N soft except one must now minimize: ( ) ( ) ( ) ∫ - ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ - - + ε κ 1 0 2 2 2 1 1 da a J a a n , (A-3) where ( ) ( ) ( ) () ∫ = π π 2
It should be mentioned that the same approach can be used when the length of the electrode is smaller than L . However, the resulting values of 0 n are no longer integers or half-integers (except when the length of the electrode goes to zero) and, as a consequence, determining the pull-in voltage or the open-loop response of the system does not boil down to finding the roots of low-order polynomials.

1

 1

  respectively the equivalent complex gains of the Duffing nonlinearity, of the electrostatic nonlinearity and of the actuation nonlinearity. Straightforward integration leads to the following expressions:

  )It is interesting to study the dependency of osc ω on osc A . It is notable that, depending on the system's parameters, there may exist an optimal value of

6 . 1 .

 61 Simulation of the 1-DOF modelThe transient simulation of (11) is conducted with Matlab/Simulink. In this environment, the pulsed force necessarily has a finite duration p τ , which introduces a small amount of phase delay into the feedback loop. To determine the influence of this non-ideality, two sets of simulations are run: one with 207 . The oscillation amplitude and the oscillation pulsation are determined at the end of the simulation,

,:

  the quadratic distance between ( ) Replacing the integrals appearing in the cost function with numerical approximations, it is possible to use the Gauss-Newton algorithm, for example, to find ε n ˆ and ε κˆ for different values of ε . The values of 0 n and 0 κ are then calculated with an extrapolation method as the limits of sequences ε n ˆ and ε κˆ. This leads to 500 . the corresponding relative error is smaller than 1% across the whole gap (figureA1) and close to 0% for 1 good results and leads to simpler analytical developments.

  a good approximation of ( ) a J

Figure 1 .

 1 Figure 1. Top-view of the mechanical structure (a) (not to scale) and diagram of the closed-loop actuation and detection scheme (b).

Figure 2 .

 2 Figure 2. Typical aspect of

Figure 3 .

 3 Figure 3. max a vs. δ , for different values of γ . The plain black curve corresponds to (20).

Figure 4 .

 4 Figure 4. Block-diagram representation of the capacitive pulse-actuated clamped-clamped beam, operated in closed-loop.

Figure 5 .

 5 Figure 5. Block-diagram representation of the capacitive pulse-actuated clamped-clamped beam, operated in open-loop.

Figure 6 .

 6 Figure 6. Amplitude (top) and phase (bottom) frequency response of the open-loop system, with 5 . 0 = γ

Figure 7 .

 7 Figure 7. SEM of one of the accelerometer structures developed in the M&NEMS project.

Figure 8 .

 8 Figure 8. osc ω vs. oscA , as predicted by (30). The relationship between osc A and f is linear.

Figure 9 .

 9 Figure 9. Relative error on the values of oscA and osc ω predicted by (30) and obtained by simulation of the 1- DOF model.

Figure 10 .

 10 Figure 10. Relative error on the values of osc A and osc ω predicted by (30) and obtained by simulation of the 3- DOF model.

Figure

  Figure A1. ( ) a I (top left), ( ) a J (top right) and relative errors (bottom) vs. amplitude.

Figure A2 .

 A2 Figure A2. Relative error between ( ) a J and ( ) a n J , , ˆ0 0 κ

  relative error is smaller than 5% across the whole gap.

	0 . One may round 0 sin 1 A I a J n to 1 and 0 2 sin dt t t . (A-4) κ to 0.182 without deteriorating the results (a 99 . 0 = 182 . relative error of less than 2% for We find 004 . 1 ˆ0 ≈ n and 0 ˆ0 ≈ κ a , figure A1). Choosing a slightly smaller value of 0 κ (say, 160 . 0 ˆ0 = κ ) yields a better approximation for moderate values of a , whereas a larger value of 0 κ ( 221 . 0 ˆ0 = κ ) improves the accuracy at very large amplitudes but reduces it for moderate a (figure A2). Any value of 0 κ in the range ] 221 . 0 , 160 . 0 [ 221 . 0 ˆ0 = κ yields the best results for 1 ≈ a can be explained analytically. Starting from ensures that the The fact that ( ) a I , 041 . 0 , 5 . 1 ˆ-

& & & &


Benoit Lange
email: benoit.lange@lirmm.fr

Nancy Rodriguez
email: nancy.rodriguez@lirmm.fr

Exploring the use of Reeb graphs in LOD generation

Keywords: I.2.10 [Vision and Scene Understanding]: Extraction of Reeb graph in 3D model Level of detail, LOD, Reeb graph, skeleton, Morse theory, OFF, QEM, Garland

Level of detail was one of the fundamental solutions of simplification 3D scene. But some issues with this solution are present. For example the topology of the mesh is not kept. In this paper we introduce the Reeb graph in level of detail. Reeb graph is a topological extractor for graph. We extract the most important point of a mesh and use them as LOD. These points are representing by a skeleton call Reeb Graph. We make a related work on LOD, skeleton extraction and extraction of Reeb graph. After, we introduce our work. We define augmented LOD or LOD+. We present the results of the different solutions that we implement like: LOD 0, LOD +, and measure between LOD. LOD 0 is the lowest LOD in a scene. LOD + consist by applying a skeleton on a 3D mesh at low resolution. And measure different between two match allow to compare different solutions of simplification.

INTRODUCTION

The modern graphics hardware is able to display thousands of polygons but they are limited in number. It is needed to simplify the 3D scene. Level of detail (LOD) is one of the most important features in graphics computation. This method is used for render large 3D model. In 1976, James Clark introduces the LOD concept [START_REF] Clark | Hierarchical geometric models for visible surface algorithms[END_REF]. It was use for produce different resolution for model. This is done by reducing the numbers of vertices and edge, this solution allow to understand easily shape of the mesh.

The most common use in LOD is global or local simplification. This kind of simplification is done on the mesh vertices or edge. But, in this paper, we use skeleton for LOD generation. Skeletons were mostly use in animation, classification, shape analysis, compression or collision detection. Skeleton can be compute or define by designer. Here, we will compute it using Reeb Graph extraction. It was introduce in 1946 by George REEB. [START_REF] Reeb | Sur les points singuliers d"une forme de Pfaff complètement intégrable ou d"une fonction numérique[END_REF] explain how to extract topology of a graph. It was an extension of Morse theory. But we can use this extractor for 3D model, it allow computing the most important point of the mesh.

This paper is decomposed in three parts.

 Related work: we introduce LOD, extraction of skeleton and Reeb graph.

 Reeb graph in LOD: we explain our four uses of Reeb graph in LOD filed.

 Results and discussion.

Related work

Many algorithms have been produce for solving the low possibility to render large meshes. David Luebke makes some survey on LOD solutions [START_REF] Luebke | A survey of polygonal simpli_cation algorithms[END_REF][START_REF] Luebke | A developer's survey of polygonal simpli_cation algorithms[END_REF]. Filter shape is the solution to simplify meshes. It"s exist three main filtering solutions.

The first solution is space division. It was use with octree, quadtree. It"s a space division of the mesh in a tree data structure. This kind of structure allow to define witch object is seen by the camera. It"s allowed to compute easily and quickly object in the scene.

Another solution is object visibility: painter algorithm or Z buffer. This two techniques use occlusion of object for hide none visible face.

The last solution study was Level of detail (LOD). It was introduced by James Clark in 1976 [START_REF] Clark | Hierarchical geometric models for visible surface algorithms[END_REF]. This solution uses many different resolutions of 3D in the scene. It necessary to define some distance between the camera and the mesh and after the algorithm display the right low resolution mesh. This part will be divided in three parts.

 Level of detail

 Skeleton extraction

 Skeleton extraction with Reeb graph

Level of detail

More and more project produce better 3D models like: the digital Michelangelo project that produce a mesh with 2 billion polygons. Actual graphics card aren"t able to display this kind of mesh. And it"s needed to simplify them. Figure 1 show different model resolution according to the distance. LOD was often release by the designer but some algorithms can produce them with keeping some important features. It exist two kinds of LOD: continues and discrete.

Discrete LOD produce the simplified model in preprocess. It was the most use solution due to the computation cost. This kind of algorithm is not in real time. We choose the number of vertices and edge and the algorithms produce models. The main issue with these algorithms is the popping effect. This effect is a visual issue: the difference between two resolution mesh is visible if the number of face is too different.

On Google earth, we can see it when we look 3D buildings. Some solutions exist to reduce this effect. We call this morphing. It consist by add vertices in real time when the resolution of the mesh change.

The other solution is continuous LOD. This method was firstly defined by Luebke. It uses a view dependent simplification to produce a mesh. This solution simplifies the entire scene instead of a individual model. These algorithms are in real time and solve the problems of popping. But the computation is expensive. And simplification in real time is hard to be done.

Then in LOD it exist two types of simplification local and global.

Local simplification is use in [START_REF] Schroeder | Decimation of triangle meshes[END_REF][START_REF] Turk | Re-tiling polygonal surfaces[END_REF][START_REF] Borrel | Multi-resolution 3d approximation for rendering complex scenes[END_REF][START_REF] Low | Model simpli_cation using vertex-clustering[END_REF][START_REF] Garland | Surface simpli_cation using quadric error metrics[END_REF][START_REF] Rau | Lod generation for 3d polyhedral building model[END_REF]. The algorithms remove, collapse, split or swap edge. These simplifications are done locally on the mesh. This can produce some hole … [START_REF] Schroeder | Decimation of triangle meshes[END_REF] is one of the first algorithm to simplify mesh. He produces multiple removing of vertex on the mesh. He done is job in three steps. Each vertex can be removing but the algorithms give some width to each one in function of their importance: complex, simple, boundary … The limitation of this simplification is the poor fidelity of the mesh. [START_REF] Turk | Re-tiling polygonal surfaces[END_REF] is a decimation solution with a sampling method.

The algorithms put on the mesh vertex, who the user defines, spring are applied on each new vertex by following the original mesh. And the method mutual tessellation produces the final model.

Many other algorithms exist and preserve topology of the mesh or some features like curves … Another way to simplify mesh is global simplification. This approach is more recent. It allows simplifying a global mesh. For example simplifications are done on the entire surface and not locally like the previous one. We can find this technique in [START_REF] Luebke | A survey of polygonal simpli_cation algorithms[END_REF][START_REF] He | Voxel based object simpli_cation[END_REF][13][14][START_REF] Eck | Multiresolution analysis of arbitrary meshes[END_REF][START_REF] Hoppe | Progressive meshes[END_REF]. Many solutions exist to simplify the global mesh: with hierarchical model, by voxel filter, by envelope simplification or base mesh.

[14] is one of this solution. It consists by reducing the number of vertices by producing two envelopes:

 Inside the mesh  Outside the mesh They are used for simplify model. The algorithm applies a simplification on the mesh. If the vertex goes out an envelope, he was bringing back between them. The big problem with this solution is the computational cost and the difficulty to extract the two envelopes. Here the method produce a base meshes and adds vertices incrementally. This can produce a better resolution for a mesh and keep the main important point of the mesh.

After this related work on LOD, we will talk about Skeleton extraction. This method is used for many problems.

Skeleton extraction

Skeleton extractions have different use graphical computation. They are use to work on the mesh with make hard computation due to the numbers of vertex.

In this part we will present the usage of skeleton extraction, some algorithms to produce this object and at last the Reeb extraction.

Usage of Skeleton extraction:

For extracting bones of a mesh, two ways are use: designer produces it by their hand or an algorithm can compute it. After the usage of this kind of structure we will present method to produce bones of the mesh.

It exist at least five usage of skeleton extraction:

 Mesh classifier  Mesh animation  Shape analysis  Compression mesh  Collision detection
The first one was classification. Here we use skeleton to recognize propriety of an object. Some algorithms produce a mesh classifier and can identify some 3D model without any human action. It allow to index mesh libraries.

A second way to use bones extraction is animation. This technique is use like bones in human body. Some part of the mesh where assign to a bone and they move like it. This solution is often create by designer because they draw there mesh around a skeleton. But some models don"t own their movement axis, and it was needed to extract them.

Mesh compression is another solution of skeleton extraction. This solution is use for transfer mesh in a little file. The skeleton was the non compressible data of the mesh and each vertex are encoded in the tree.

Extracting skeleton can help for collision. Often, we use hull, box or sphere to detect the collision. Algorithm like Qhull use voronoï diagram to produce the hull. The problem with this solution is the non realistic collision issue. An animal cannot be represent by a box due to is tree structure. So use a skeleton is more realistic collision detection. The width of the bone is a problem due to the thickness of it.

Finally, the last use is shape analyze. This method allows extracting some important feature on the mesh: we can identify the curve or something else. The main use for extract them is REEB algorithm. In this paper we will talk about this solution.

In the next part, we will talk about different way to extract skeleton.

Skeleton extraction:

Many solutions exist for extract the mesh skeleton.

Here we will talk about two methods.

[17] is middle axe geodesic extraction. He produces the middle axis of each geodesic surface of the mesh. After this we compute the middle axis by align each geodesic axe. This method is one of the first who the mesh position is not important. With this solution we can identify some shape of the mesh like hand.

[18] is another solution who compute the middle axis of the mesh by using spring on each vertex. After a few iterations, this algorithm aligns bones and aligns important axis like shoulders or something else. This kind of algorithms is mostly use for animals.

Other solution for extracting skeleton is topological extraction. For doing this we use [START_REF] Clark | Hierarchical geometric models for visible surface algorithms[END_REF] and we extract le main important points. This graph is call Reeb graph.

He was a evolution of Morse theory on manifold graph. This solution extract singular point of a shape .

Morse theory was introduced in 1925 by Marton

Morse. This theory infers manifold"s topological invariant from a Morse function. In [START_REF] Tierny | Analyse topologique et géométrique de maillages 3D pour l'extraction de squelette[END_REF], Tierny give an mathematical translation of the Reeb graph theory.

In fact Reeb graph are compose from nodes that are the critical point of the structures. And they are the nodes of the reeb graph. Critical points are extracting from a special function like diameter of a graph. For each point of the long way of the graph, we compute critical points by using gradient function.

It exist three critical points to identify. The first one was bijection point. It"s a point who the next step, on the diameter way, separate himself to produce a branch. Another important point is merging point. The branch of the tree becomes a single point. And finally the next step try to identify the end point of the mesh.

Reeb graph in LOD:

In this paper we will introduce four new methods who use Reeb Graph: LOD 0, LOD +, Drive LOD and measure of LOD.

We use harmonic extraction of Reeb graph but these techniques have some problems.

In fact mesh need to be manifold and it"s not easy to find this kind of models.

Another issue was file format. The simplification use for simplify is QEM from Garland who use SMF files and the harmonic extraction of skeleton use OFF files. We need to use a file convertor to produce correct files.

Skeleton have to type define, so we need to create a special file format for it. We define a format like OFF file call Skel file.

SKEL V 0 0 0 0 V 1 1 1 1 V 2 0 1 1 E 0 1 E 1 2
The V lines are for vertices and the E line are for edges.

Now, we will talk about the first part of our contribution: LOD 0. All of our new usage of Reeb graph uses several meshes: a hand, a boy, a chair, a duck, a horse and a CAD object. Before compute our algorithms, we can identify 2 kind of objects: with many important topological point (ex: hand, boy, chair or horse) and the other one (ex: duck, CAD).

The first class can be assimilating like a tree.

LOD 0:

Level of detail has a big issue. If the algorithms make too hard simplification the shape wasn"t keep, and the topology of the mesh become totally wrong. The model become impossible to understand and is usage is totally useless. Our Solution consists by using a very simple mesh for view the model at a far distance. This method was call LOD 0 or LOD max. The first are good on visualization but the width of the hand is not conserved. Two problems are popping from this solution.

The first one was width. Figure 7 show that the width of the edge is important. So we decide to use some tricks to dress the tree. We use three "clothes": box, cylinder-ball, and cylinder. Cylinder and balls produce good result in shape recognition but are cost expensive. Other solutions are not efficient: they don"t respect the shape of the mesh and produce some artifact like holes.

The second problem was on the kind of mesh. For memory we define two mesh classes: tree and none tree. The trees meshes have really good results, we can find 5 fingers on the mesh while QEM merge them. But with none tree meshes, we cannot imagine the shape. Figure 8 The skeleton of the CAD piece is representative of the mesh but without "clothes" the shape of the mesh is not identifiable.

In the next part we will talk about the next solution of our solution.

LOD +:

A new feature in LOD domain is Augmented LOD or LOD +. He was the same solution like augmented reality but we set the skeleton on the simplified mesh. This solution allows keeping some important features: for example when we compute the hand with QEM, fingers disappear or merge. With this solution we don"t lost anymore the important stuff of the mesh. But two problems income:

 The first was due to the width of the mesh like the first solution.  The second issue was on none tree meshes.

We don"t find the important features of the mesh and this technique was not useful. The next step will introduce quickly a solution that we haven"t implement.

Assisted LOD:

This solution consist to drive algorithm in chooses of simplification.

Many algorithms use their own solution for drive the simplification. We imagine some solution on many algorithms.

For like [START_REF] Luebke | A survey of polygonal simpli_cation algorithms[END_REF]: the clustering grid can be compute by using Reeb graph. This solution allows the simplification to be done on the important features of the mesh. Topology of the mesh was conserved and the simplified mesh will have a good visual aspect.

In [START_REF] Garland | Surface simpli_cation using quadric error metrics[END_REF], we can compute the graph each time that we remove an edge. This allows measuring the topological modification of the simplified mesh at each modification. One of the main issues of this solution is the computation cost. Reeb graph is O(N(log(N)), if we apply the same algorithms at each topology modification, the cost will be too expensive. Some heuristics can be used to not recompute the entire graph. But many other algorithms, who didn"t talk, can use Reeb graph. The last exploration was measure between two mesh.

Measure of a mesh and his simplification:

In this part, we will talk about a new solution to measure different between an original mesh and is simplification. Two previous solutions exist.

The existent solution use simplification method or sampling techniques. Here we want to compute the Reeb graph of a mesh and is simplification. And we compare the two graphs. More graphs are similar more the simplification have been efficient. The first step of algorithm is to compute the two skeletons. After, we try to approximate a theoretical Reeb graph with some parameters: numbers of vertices, location of vertices. We compute the new graph and compare to the simplified tree. With this solution we can identify witch cluster contain the most problems. The other solutions of measure don"t identify the branches that have troubles.

Other ways to identify problems are possible but we don"t study them in this paper.

Results and discussion

We introduce four approaches to use Reeb graph in LOD. We have implemented three of them. And results are quite good.

We find that some meshes aren"t usable with some techniques: none tree mesh. But other models have many benefits of our solutions.

For LOD 0, the solution is good for tree meshes but its need to find a way to wear the skeleton because it"s too thin for show at far distance.

For LOD +, we can restore some important part of the meshes. Little details that disappear with simplification method are keeping. And finally, measure LOD is efficient. She allows matching the branch with the biggest issues. Our technique is not the best but she gives some interesting results.

On our future work, we want to study new implementation of Reeb graph [START_REF] Tierny | Graphes de Reeb de Haut Niveau de Maillages Polygonaux 3D[END_REF].

We want to implement the assisted LOD witch many algorithms like: [START_REF] Luebke | A survey of polygonal simpli_cation algorithms[END_REF][START_REF] He | Voxel based object simpli_cation[END_REF]13]. This allows to computer better low resolution mesh instead of using their metrics. The topological cost can be the key of mesh simplification to produce better LOD models.

Another solution is to implement other solution to compute metrics between a object and the simplified mesh.

Figure 1 .

 1 Figure 1. Different resolution of a mesh

Figure 2 Figure 3 .

 23 Figure 2. Original Mesh

Figure 4 .

 4 Figure 4. Simplification envelope

Figure 5 .

 5 Figure 5. Reeb graph of a simple shape

Figure 6 .

 6 Figure 6. Hand in full resolution

Figure 6

 6 Figure 6 is the original mesh and figure 7 present is Reeb graph in the scene.

Figure 7 .

 7 Figure 7. Hand Reeb graph at far distance

  and 9 show the CAD piece and the Reeb graph of it.

Figure 8 .Figure 9 .

 89 Figure 8. CAD piece at full resolution

Figure 10 .

 10 Figure 10. Hand at low resolution Figure 11 present the LOD + who give back to the mesh the fingers.

Figure 11 .

 11 Figure 11. LOD +

ACKNOWLEDGMENTS

We thank every people who help use during the study. And we want to thanks Urbsim to give the idea of use Reeb graph in LOD.