T G Ritto
email: tritto@mecanica.ufrj.br

C Soize

F A Rochinha

Rubens Sampaio
email: rsampaio@puc-rio.br

Dynamic stability of a pipe conveying fluid with an uncertain computational model

Keywords: fluid-structure interaction, uncertainty quantification, stochastic dynamics, stochastic stability analysis, reliability analysis

) to analyze the stability and reliability of the stochastic system. The Euler-Bernoulli beam model is used to model the pipe and the plug flow model is used to take into account the internal flow in the pipe. The resulting differential equation is discretized by means of the finite element method and a reduced-order model is constructed from some eigenmodes of the beam. A probabilistic approach is used to model uncertainties in the fluid-structure interaction. The proposed strategy takes into account global uncertainties related to the noninertial coupled fluid forces (related to damping and stiffness). The resulting random eigenvalue problem is used to analyze flutter and divergence unstable modes of the system for different values of the dimensionless flow speed. The numerical results show the random response of the system for different levels of uncertainty, and the reliability of the system for different dimensionless speeds and levels of uncertainty.

Introduction

This paper extends the deterministic stability analysis proposed by [START_REF] Paidoussis | Dynamic stability of pipes conveying fluid[END_REF] of a pipe conveying fluid. The present work deals with a probabilistic model that takes into account uncertainties induced by modeling errors that arise due to physical simplification introduced in the deterministic model, as it will be explained latter. Slender flexible tubes with internal flow or pipes conveying fluids are present in a number of applications, such as microsystems, biological devices, drill-strings and heat exchangers. See, for instance, [START_REF] Ritto | Nonlinear dynamics of a drill-string with uncertain model of the bit-rock interaction[END_REF] for nonlinear dynamics of a drill-string, [START_REF] Rinaldi | Dynamics of microscale pipes containing internal fluid flow: Damping, frequency shift, and stability[END_REF] for microscale resonators, [START_REF] Soltani | Vibration and instability of a viscous-fluid-conveying single-walled carbon nanotube embedded in a visco-elastic medium[END_REF] for nanotubes with viscous fluid/ Kelvin-Voigt model, and [START_REF] Gosselin | Dynamical stability analysis of a hose to the sky[END_REF] for the dynamic stability of a hose to the sky.

Typically, the standard configuration is a straight tube mounted over supports carrying a steady flow with a constant velocity. If the velocity is kept below a certain threshold, ambient perturbations might entail low amplitude vibrations of the structure around a standard configuration. Above the critical speed, the coupled system, constituted by the structure and the flow, might undergo large vibrations and complex nonlinear dynamical responses. The understanding of such an unstable behavior is required for improving the design of the system and mitigate damage effects. In addition, the large diversity of the dynamical response renders this problem quite attractive for theoretical and numerical studies.

Computer models are nowadays widely used in the design and analysis of standard engineering systems. Many critical decisions are taken based on computational simulations. Despite the consolidation of powerful and reliable methods leading to small numerical errors originated by discretization techniques, the extension of this common practice to more critical systems is hindered by the presence of inevitable uncertainties associated with the modeling. Statistical fluctuations around nominal values of parameters, nonidealized initial and boundary conditions, or production tolerances might entail a large variation on the output of the simulations. So, in order to improve the reliability of predictions, those uncertainties must be taken into account. Here, the consideration of the uncertainties is carried out within a probabilistic framework, hence, design criteria are based on failure probabilities and reliability analysis, such as done by [START_REF] Ghanem | Spectral stochastic finite-element formulation for reliability analysis[END_REF] and [START_REF] Schueller | A critical appraisal of reliability estimation procedures for high dimensions[END_REF]. The problem of a pipe with internal flow can be modeled using a tridimensional nonlinear elastic model for the structure together with Navier-Stokes equations for the fluid; see [START_REF] Bathe | Finite element analysis of fluid flows fully coupled with structural interactions[END_REF] and [START_REF] Bathe | Finite element developments for general fluid flows with structural interactions[END_REF]. But, usually, simplified physics is introduced. In Piet-Lahanier and Ohayon (1990), for instance, a beam model is used for the pipe, and a compressible and viscous fluid is considered. The use of simplified models, such as the one proposed by [START_REF] Paidoussis | Dynamic stability of pipes conveying fluid[END_REF] and reused in this paper, makes feasible the analysis of a significant number of scenarios. However, it is clear that the use of kinematic reductions introduces modeling errors.

The stability of dynamic systems has been extensively studied from a deterministic perspective; see, for instance, [START_REF] Guckenheimer | Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields[END_REF] and [START_REF] Nayfeh | Applied Nonlinear Dynamics[END_REF]. There are also some works in the literature dealing with the stochastic stability analysis. For example, in [START_REF] Lin | Probabilistic Structural Dynamics[END_REF], the concepts of almost-sure stability, stability in probability and stability in the mth moment are defined and used for analyzing stochastic dynamical systems submitted to time-varying loads.

The stability of pipes with internal flow was deeply investigated by [START_REF] Paidoussis | Fluid-Structure Interactions: Slender Structures and Axial Flow[END_REF]. In [START_REF] Paidoussis | Dynamic stability of pipes conveying fluid[END_REF], some historical review of the subject up to the time the article was written is done, and the instabilities due to divergence and flutter phenomena are discussed for steady-state flows and harmonically perturbed flows. This was the first time that the coupledmode flutter (simply called flutter) was noticed for this kind of coupled system. In [START_REF] Ariaratnam | Dynamic stability of pipes conveying pulsating fluid[END_REF], an analytical method (devoted to the stability of pipes with perturbed internal flow) is proposed by using the method of averaging and the Floquet-Lyapunov theory. In [START_REF] Ibrahim | Overview of mechanics of pipes conveying fluids-part I: Fundamental studies[END_REF], a more recent review of the mechanics of pipes conveying fluid can be found in which more than four hundred references are given.

There are few investigations related to uncertainty for stability of a pipe with internal flow. In Ariaratnam and Namachchivaya (1996b), random velocity fluctuations were considered; the authors used the averaging method and the Floquet theory. In [START_REF] Ganesan | Vibration and stability of fluid conveying pipes with stochastic parameters[END_REF], system-parameter uncertainties are taken into account. Some statistics of the critical flow ve-locity are analyzed for a stochastic modeling of the elasticity modulus and of the mass per unit length of the structure. In [START_REF] Yigit | Active control of flow-induced vibrations via feedback decoupling[END_REF], the flutter stability of a cantilever pipe conveying fluid is considered, where active control is used to suppress the structural vibration. The present paper aims to analyze the dynamic stability of a pipe conveying fluid including uncertainty in the fluid-structure interaction model, as explained in the sequence.

The fluid-structure interaction phenomena are responsible for the existence of unstable modes of the dynamical system and, therefore, play a central role in the present analysis. In order to endow the model with an improved capacity of describing the fluid-structure coupling, a probabilistic approach which has the capability to take into account modeling errors should be used. In this paper, we propose to use the nonparametric probabilistic approach introduced by [START_REF] Soize | A nonparametric model of random uncertainties for reduced matrix models in structural dynamics[END_REF][START_REF] Soize | Stochastic Models of Uncertainties in Computational Mechanics[END_REF], which is a method to take into account model uncertainties induced by modeling errors.

In computational dynamics and in computational fluid-structure interaction such as aeroelasticity and vibroacoustics, the nonparametric probabilistic approach is an alternative method to the output-prediction-error method which allows modeling errors to be taken into account at the operators level of the computational model by introducing random operators and not at the model output level by introducing an additive noise. It should be noted that such an approach allows a prior probability model of uncertainties to be constructed even if no experimental data are available. The nonparametric probabilistic approach is based on the use of a reduced-order model and the random matrix theory. It consists in directly constructing the stochastic modeling of the operators of the mean computational model. The random matrix theory and its developments in the context of dynamics, vibration and acoustics is used to construct the prior probability distribution of the random matrices modeling the uncertain operators of the mean computational model. This prior probability distribution is constructed by using the Maximum Entropy Principle introduced by [START_REF] Shannon | A mathematical theory of communication[END_REF]; [START_REF] Jaynes | Information theory and statistical mechanics[END_REF][START_REF] Jaynes | Probability Theory: The Logic of Science[END_REF], in the context of Information Theory, for which the constraints are defined by the available information. Section 4 is devoted to the nonparametric probabilistic model, where the statistically dependent random matrices which are derived from the fluid coupling model are generated thanks to the introduction of an adapted scheme. Consequently, with the aid of the nonparametric probabilistic approach, the uncertainties are globally modeled and modeling errors are taken into account. The nonparametric probabilistic approach of uncertainties has been applied in numerous different areas to cope with uncertainties arising from different sources along the modeling process; see, for instance, [START_REF] Chen | Probabilistic approach for model and data uncertainties and its experimental identification in structural dynamics: Case of composite sandwich panels[END_REF]; [START_REF] Durand | Structural-acoustic modeling of automotive vehicles in presence of uncertainties and experimental identification and validation[END_REF][START_REF] Ritto | Nonlinear dynamics of a drill-string with uncertain model of the bit-rock interaction[END_REF]. This probabilistic approach was also used by Ritto et al. (2011) to model uncertainties in the stiffness of a piping conveying fluid but no stability analysis was done. Since the present work has its particularities, a novel idea is proposed to adapt the procedure of the nonparametric probabilistic approach.

Considering the randomness of the computational model, the stability analysis relies on solving a stochastic eigenvalue problem, [START_REF] Croquet | Étude des dispersions et incertitudes en optimisation et dans l'analyse des valeurs propres[END_REF]. A number of results are presented with emphasis on the statistical characterization of the eigenvalues. (1) A sensitivity analysis is carried out as follows. For different values of the dimensionless fluid speed and for different levels of uncertainty, the confidence regions of the eigenvalues and of the stochastic response are constructed using the Monte Carlo method, [START_REF] Rubinstein | Simulation and the Monte Carlo Method[END_REF], as stochastic solver. (2) A reliability analysis is then performed. The probability of being in the divergence region or in the flutter instability region is computed. In addition, a control surface is drawn showing, for different values of the dimensionless fluid speed and of the levels of uncertainty (u,δ), how severe is the flutter instability. This paper is organized as follows. The deterministic system is depicted in Sections 2 and 3, where the equations, the reduced-order model and the stability analysis are presented. The probabilistic model is developed in Section 4. Finally, the numerical results are analyzed in Section 5 and the concluding remarks are made in Section 6.

Deterministic model

Although the deterministic formulation of a pipe conveying fluid is becoming relatively classical, [START_REF] Paidoussis | Fluid-Structure Interactions: Slender Structures and Axial Flow[END_REF], for the sake of completeness, some details are presented in this section. Consistently with the aim of studying the stability of the coupled system, the rectilinear configuration of the coupled system is sketched in Fig. 1. The Euler-Bernoulli beam theory for an undamped simply supported beam is adopted. Therefore, the partial differential equation governing the dynamics of the structure is written as; see, for instance, [START_REF] Inman | Engineering Vibration[END_REF]:

m ∂ 2 v(x, t) ∂t 2 + EI ∂ 4 v(x, t) ∂x 4 = f (x, t) , (1)
with appropriate boundary and initial conditions, where v is the transversal displacement, L is the length of the beam, m is the mass per unit length, E is the elasticity modulus, I is the area moment of inertia, f is the external force per unit length, x and t are the position and time parameters that are defined in the intervals [0, L] and [0, T ]. To model the internal flow, the plug flow model, [START_REF] Paidoussis | Fluid-Structure Interactions: Slender Structures and Axial Flow[END_REF], is used. Let r = xi + yj be a point measured from the origin, where (i, j)is fixed in an inertial reference frame, with j the unit vector in the transverse direction and where i is the unit vector perpendicular to it, aligned with the axial direction. Due to the assumption that the beam undergoes only small displacements, the velocity of a particle of the fluid can approximately be expressed as

v f ≃ Ui + U ∂v ∂x j + ∂v ∂t j , ( 2
)
where U is the speed of the fluid with respect to the structure and where the rotation, ∂v/∂x,of the beam cross sectional area is assumed to be small. Hence, the fluid acceleration is given by

a f = Dv f Dt = ∂v f ∂t + (v f •∇)v f . (3)
Considering a constant speed U (i.e. dU/dt=0) yields

a f = ∂ 2 v ∂t 2 + 2U ∂ 2 v ∂x∂t + U 2 ∂ 2 v ∂x 2 j . (4)
The force per unit length is obtained multiplying the acceleration by the fluid mass per unit length M f ,

f f = M f ∂ 2 v ∂t 2 + 2M f U ∂ 2 v ∂x∂t + M f U 2 ∂ 2 v ∂x 2 j . (5)
Substituting this force in Eq. ( 1) yields,

(m + M f ) ∂ 2 v ∂t 2 + 2M f U ∂ 2 v ∂x∂t + M f U 2 ∂ 2 v ∂x 2 + EI ∂ 4 v ∂x 4 = f . (6)
The above equation shows that the fluid contributes to the mass, damping and stiffness of the system. Equation ( 6) is discretized by means of the finite element method. Introducing the local coordinate ξ = x/l e of a finite element, E e , in which l e is the element length, for x in E e , the finite element approximation of the displacement field v(x, t) is written as v (e) (ξ, t) = N(ξ)u (e) (t) in which N is the shape Hermitian function and where u (e) = [v 1 ∂v 1 /∂ξ v 2 ∂v 2 /∂ξ]i s the element displacement vector. The corresponding mass, fluid damping, fluid stiffness and bending stiffness element matrices are then given by M (e) = (m

+ M f ) 1 0 N T Nl e dξ, C f (e) = 2M f U 1 0 N T N ′ dξ, K f (e) = -M f U 2 1 0 N ′T N ′ 1 l e dξ and K b (e) = EI 1 0 N ′′T N ′′ 1 l 3 e dξ.
The assembly of the finite element matrices yields the following computational model,

M ü(t) + C f (U) u(t) + (K b + K f (U)) u(t) = f(t) , ( 7
)
where M is the positive-definite (m × m) mass matrix, K b is the positive-definite (m × m) bending matrix and K f (U) is the negative-definite (m × m) fluid stiffness matrix (for U 0), u(t)i st h e response vector and f(t) is the load vector. In general, the (m × m) real matrix C f (U) is not symmetric and, hence, can be written as the sum of a symmetric positive-definite matrix with a skew-symmetric one. Due to the boundary conditions of the present application, which is simply supported, C f (U) is skew-symmetric. However, we still call it damping matrix, even though it does not entail energy dissipation. In the frequency domain, Eq. ( 7) can be written as

(-ω 2 M + iω C f (U) + K b + K f (U)) û(ω) = f(ω) , (8)
where û(ω) is the m-complex response vector, f(ω)isthem-complex load vector, ω is the circular frequency and i = √ -1. It should be noted that the application f → û(U,ω) is linear, but the mapping U → û(U,ω) is nonlinear. Indeed, for the main goal here consisting on carrying out the stability analysis of the coupled system based on the numerical model, the mapping of interest is implicitly defined by Eq. ( 8) relating the speed U to the eigenvalues, that is U → λ i (U) (i = 1, 2,..).

Reduced-order computational model and stability analysis

A reduced-order computational model, [START_REF] Holmes | Turbulence, coherent structures, dynamical systems and symmetry[END_REF]; [START_REF] Kerschen | The method of proper orthogonal decomposition for dynamical characterization and order reduction of mechanical systems: An overview[END_REF]; [START_REF] Murthy | Uncertainty-based experimental validation of nonlinear reduced order models[END_REF], is going to be constructed for the system given by Eq.( 8). In structural dynamics and coupled systems, the construction of a reduced-order computational model is very useful (1) to get a model adapted to the prediction of the physical phenomena under consideration, (2) to facilitate the analysis of the predictions (in particular for stability analysis), (3) to reduce the numerical cost (in particular, in the context of stochastic dynamics for which a Monte Carlo method is used as stochastic solver) and (4) to make feasible the nonparametric probabilistic approach explained in details in the next Section.

Since the dynamical system is linear, the family of eigenmodes is generally considered the most efficient basis for the construction of a reduced-order computational model. When the dynamical system is nonlinear, a reduced-order model can always be constructed using an adapted basis. The efficiency of such a basis will then be quantified by the numerical cost of its construction and by the speed of convergence with respect to the dimension of the reduced-order computational model. In this context, the family of eigenmodes related to a linear dynamical system associated with the nonlinear one can still be efficient. In that way, the projection basis is independent of the external forces applied to the dynamical system. In other cases, it can be interested to use proper orthogonal decomposition (POD) to construct the basis, but such a method requires to prior solve the problem in the physical coordinates and, in addition, the constructed basis depends on the external forces applied to the dynamical system. Introducing the homogeneous equation associated with Eq.( 8):

(-ω 2 M + iω C f (U) + K f (U) + K b ) û(ω) = 0 . ( 9
)
The following generalized eigenvalue problem, corresponding to Eq. ( 9) for U = 0, is introduced to construct the reduction basis:

(-ω 2 M + K b ) φ = 0 . ( 10
)
For n < m, let Φ = [φ 1 φ 2 ... φ n ] be the (m × n) real matrix composed by the first n eigenmodes. With this choice, (1) we guarantee an appropriate basis (M and K b are positive-definite matrices) and (2) it is not necessary to recompute the basis for each value of U. The displacement is then written as û(ω)≃Φ q(ω) and the reduced-order computational model is given by

(-ω 2 M r + iω C r (U) + K r (U)) q(ω) = fr (ω) , (11)
where fr (ω) = Φ T f(ω) and the reduced matrices are given by

M r = Φ T MΦ , K r (U) = K fr (U) + K br , K fr (U) = Φ T K f (U)Φ , K br = Φ T K b Φ , C r (U) = Φ T C f (U)Φ , (12)
The normalization of the eigenmodes is chosen such that M r = I (where I is the (n × n) identity matrix) and, consequently,

[Φ T K b Φ] ij = δ ij ω 2 i ,
in which δ ij is the Kronecker symbol. The stability of the system is analyzed solving an eigenvalue problem constructed as follows. The size of the system is doubled (state space in the frequency domain) and Eq. ( 11) is rewritten as

0I -K r (U) -C r (U) B r (U) q iω q y = iω I0 0M r D r q iω q (13)
Since M r = I and introducing λ = iω, the eigenvalue problem associated with Eq.( 13) assumes the conventional form:

B r (U) y = λ y , ( 14
)
where y is the eigenmode associated with the eigenvalue λ = Re(λ) + i Im(λ) in which Re(•) and Im(•) refer respectively to real and imaginary components. If Re(λ) > 0, then the system is unstable. If Re(λ) = 0, then the system is unstable if the corresponding eigenvalues are not simple roots of the characteristic polynomial.

Probabilistic model of uncertainties

As mentioned in Section 1, the nonparametric probabilistic approach, [START_REF] Soize | A nonparametric model of random uncertainties for reduced matrix models in structural dynamics[END_REF][START_REF] Soize | Random matrix theory for modeling uncertainties in computational mechanics[END_REF], is used to take into account model uncertainties induced by modeling errors in this fluid-structure interaction problem. This probabilistic strategy is chosen because it is a way to cope with model uncertainties introduced by the different sources like, for instance, simplifications on the kinematics used in the present model; which reduce its predictive ability. The nonparametric probabilistic approach couples eigenmodes coordinates, which seems well suited to this problem, since coupling is involved. In this coupled problem, the sources of uncertainties are the following: structural uncertainties (use of Euler-Bernoulli beam theory, boundary conditions, material properties) and fluid-structure coupling uncertainties (velocity field approximation, fluid properties). In the present paper, only fluid-structure coupling uncertainties are the subject of analysis. Therefore, uncertainties related specifically to the structure and uncertainties in the mass properties or external forces are not taken into account.

In this context, uncertainties in the noninertial fluid forces, which are related to the fluid damping and stiffness matrices, are globally taken into account. Since these forces come from the same coupled model, they are statistically dependent from each other, therefore, we propose a probabilistic model that considers this fact and, more importantly find a way to generate their joint distribution.

Let us first present a quick review of the concepts supporting nonparametric probabilistic approach. Let matrix D be an operator of the dynamical problem under consideration. Matrix D is assumed belonging to a set M of (n × n) real matrices. For instance, D could be the reduced stiffness matrix K br which belongs to the set M = M + of all the positive-definite (n × n) real matrices. Assuming that the matrix D depends on a vector-valued parameter, p, and calling it by D(p). If p is uncertain, in the context of the usual parametric probabilistic approach of system-parameter uncertainties, it is modeled as a random vector p for which its probability distribution has to be constructed and which is defined on the admissible set S of parameter p. In such a case, the probability distribution of the random matrix D par = D(p) is induced by the probability distribution of p. The family of deterministic matrices {D(p), p ∈S }spans a subset M 0 of set M. Consequently, the probability distribution of random matrix D par is defined on subset M 0 and not on the entire set M. This means that the parametric probabilistic approach cannot reach any matrix in set M but only those matrices that are in subset M 0 . If there are no modeling errors, but only system-parameter uncertainties, then the "experimental" matrix, D exp , corresponding to D par , belongs to set M 0 , and the parametric probabilistic approach will be efficient because the probability distribution of random matrix D par is defined on M 0 . If there are modeling errors, D exp is outside set M 0 and consequently, the parametric probabilistic approach has not the capability to reduce the distance between D par and D exp . In order to take into account modeling errors, the nonparametric probabilistic approach consists in directly constructing the probability distribution of the random matrix D which models D on the set M and consequently, such a model allows D exp to be reached by D. Figure 2 shows the Venn diagram that summarizes this concept. Another key point is the construction of the probability distribution, also called the stochastic structure, of random matrix D such that the properties of matrix D are respected.

The methodology of the nonparametric probabilistic approach consists in (i) constructing a reduced-order model deduced from the nominal computational model (also called the mean computational model), (ii) replacing the deterministic matrices of this reduced-order model by random matrices which allows the model uncertainties to be taken into account (this is the main concept), (iii) introducing the sets of random matrices consistent with the mathematical properties of these matrices and including all the available information and (iv) constructing the probability distribution of the random matrices belonging to each set of random matrices using the Maximum Entropy Principle.

For instance, if D is positive definite, it can be decomposed as D = L T IL , where I is the identity matrix and L is an upper triangular matrix. Then, the random matrix D is constructed as D = L T GL, with G being a random matrix with mean equals to I. The stochastic structure of random matrix D is then defined by this construction and by the additional available information relative to G which can be summarized as follows for a large class of problems:

• G is a positive-definite random matrix,

•E { G} = I , •E { | | G -1 || 2 F } = c 1 , |c 1 | < +∞ , where ||A|| F = (trace{A T A}) 1/2
is the Frobenius norm and E{•} denotes the mathematical expectation. The first information says that random matrix G (and, consequently D) has values in the set M + , which is an important property to be respected. The second condition simply says that the mean value of random matrix D is equal to the nominal value, D. The last condition enforces the variance of the response to be finite, which must be true to ensure that the response has a physical meaning. The probability distribution of G is then constructed using the Maximum Entropy Principle under the constraints defined by the above available information. Its expression and the corresponding generator of independent realizations can be found in [START_REF] Soize | A nonparametric model of random uncertainties for reduced matrix models in structural dynamics[END_REF][START_REF] Soize | Random matrix theory for modeling uncertainties in computational mechanics[END_REF]. This set of random matrices G is denoted by SG + . It should be noted that other sets of random matrices, including additional available information, has also been constructed for structural dynamics, Mignolet and Soize (2008b), and that a combination of the parametric and the nonparametric probabilistic approaches (generalized probabilistic approach) is possible, [START_REF] Soize | Generalized probabilistic approach of uncertainties in computational dynamics using random matrices and polynomial chaos decompositions[END_REF].

Going back to the problem analyzed in the present paper, a stochastic structure is proposed for the noninertial forces related to the action of flowing fluid over the pipe in order to take into account model uncertainties. It should be noted that this interaction encompasses two coupling mechanisms (stiffness and damping) which have to be taken into account simultaneously, as both are governed by the flow speed. In addition, the algebraic properties of both matrices implied from the physics of the problem must be respected, and, in the present model, the deterministic matrices are not frequency dependent. Therefore, we would like to enforce that the stochastic models of these matrices should not depend on the frequency. Since we do not have enough information to characterize the statistical dependence between the two matrices in line with their strong dependency on U, we will use a single random germ for the two matrices, as will be depicted in the sequence. For these conditions, we can gather the information of both damping and stiffness matrices originated by the fluid motion in a single one A(U) as follows

A(U) = iω r C r (U) + K fr (U) , (15)
where ω r is a reference frequency. Then, A(U) is modeled by a random matrix A(U), with construction carried out as explained in {Soize05 and [START_REF] Mbaye | Robust analysis of design in vibration of turbomachines[END_REF]. The deterministic matrix A(U), which is invertible for all U 0 (since K fr is symmetric negative definite), is first decomposed by the polar decomposition, Golub ad Van Loan (1996),

A(U) = Q(U) P(U) . ( 16
)
From now on, the dependence on U will be omitted to simplify the notation. In Eq. ( 16), Q is a unitary matrix (Q * Q = QQ * = I), in which * denotes the conjugate transpose, and P is a Hermitian positive-definite matrix. This decomposition will be computed using the singular value decomposition, Golub ad Van Loan (1996): A = USV * , with P = VSV * and Q = UV * . The Cholesky decomposition can be applied to the Hermitian positive-definite matrix P = L * L.

The construction of the random matrix A consists in writing

A = QL * G L , ( 17
)
where G belongs to the set SG + of the random matrices introduced above. For a given realization of G, it is straightforward to compute the corresponding realizations of A. The level of uncertainty related to the fluid-structure interaction model is controlled by the dispersion parameter δ of matrix G, which is defined as

δ = 1 n E{||G -I]|| 2 F } 1 2 . ( 18
)
It should be noted that only one parameter, δ, controls the dispersion of the stochastic model, which considerably simplifies the identification procedure when experimental data are available, solving a statistical inverse problem, [START_REF] Soize | Probabilistic model identification of uncertainties in computational models for dynamical systems and experimental validation[END_REF]; [START_REF] Soize | Stochastic Models of Uncertainties in Computational Mechanics[END_REF]. With the stochastic structure defined by Eq.( 17), the uncertainties are globally taken into account in random matrix A and uncertainties mainly connected to the flow (kinematic simplifications, non-uniformity,...).

The only requirements are that (1) the mean value of A is equal to matrix A of the nominal model, (2) that random matrix G belongs to SG + which means that G is a positive-definite random matrix and that the random response of the stochastic reduced-order computational model will have a finite variance. The random eigenvalue problem related to Eq.( 13) is written as:

0I -K r (U) -C r (U) B r (U) y = i Ω I0 0M r D r y ( 19
)
in which K r (U) = Re{A} + K br and C r (U) = Im{A}/ω r , for ω r > 0. The reference frequency ω r could be, for instance, the mean value in the frequency band. Introducing Λ=i Ω, Eq. ( 19) can then be rewritten as

B r (U) y =Λy , ( 20
)
where the complex-valued random eigenvalue is written as Λ=Re(Λ) + iIm(Λ). Finally, the stochastic reduced-order computational model associated with Eq. ( 11) is written as

û(ω) = Φ q(ω) , ( 21
) (-ω 2 M r + iω C r (U) + K r (U)) q(ω) = fr (ω) , (22)
in which û(ω) is the random response associated with the deterministic response û(ω).

Numerical results

The structure is discretized with 40 finite elements (choice made after a convergence check) and the first four eigenmodes are considered in the analysis. To be able to do a more general analysis, some dimensionless quantities are introduced:

ζ = x/L, η = v/L, β = M f /(m + M f ), u = UL(M f /(EI)) 1/2 , f = fL 3 /(EI), ̟ = ωL 2 ((m + M f )/(EI)) 1/2 ,
where ζ is the dimensionless length, η is the dimensionless transverse displacement, β is the ratio of fluid mass, u is the dimensionless speed, f is the dimensionless force, and ̟ is the dimensionless frequency. The dimensionless eigenvalues are written as λ nd and Λ nd . There are two important dimensionless parameters which are β (mass density relation) and u (dimensionless flow speed). In the stability analysis, we set β = 0.24 and the value of the dimensionless speed u varies. The frequency band of analysis is ̟ ∈ [0, 35]. For the frequency domain analysis, the value of the dimensionless force applied to each degree of freedom of the beam is 1 × 10 -4 and the response displayed is the absolute value of the dimensionless displacement in the middle of the beam.

Deterministic analysis

This section deals with the deterministic analysis of the system, that will be compared with the stochastic analysis of the next section. The results presented are in accordance with the ones found in [START_REF] Paidoussis | Fluid-Structure Interactions: Slender Structures and Axial Flow[END_REF].

Figure 3 shows the stability charts using the real component (Fig. 3(a)) and the imaginary component (Fig. 3(b)) of the eigenvalues. As u increases, the system becomes more unstable. The first divergence mode occurs at u = π and the second one at u = 2π. The first flutter instability occurs around u = 6.29 and the third divergence mode occurs around u = 3π,a s shown in the Figure. Figure 3(b) shows that, as u increases, the system becomes less stiff, which is observed by the decreasing of the modulus of the imaginary parts of the eigenvalues considered in the analysis. Figure 4 shows the amplitude of the displacement in the frequency domain for u = {1.0, 2.5}. The response corresponding to the higher value of u is displaced to the left (system less stiff).

The results presented in this section brings no new information with respect to what has been published in the literature. It simply reproduces known behaviors of a pipe conveying fluid, [START_REF] Paidoussis | Fluid-Structure Interactions: Slender Structures and Axial Flow[END_REF], but it is necessary to produce these results to validate the deterministic computational model and to use it as a reference for the stochastic analysis. The next section analyzes how these results are modified when the uncertainties are taken into account.

Uncertainty quantification

The Monte Carlo method is used as the solver of the resulting stochastic system, and 100, 000 deterministic simulations were performed, such that the mean-square convergence of the eigenvalues was achieved. Section 5.2.1 shows the confidence regions of the analysis and Section 5.2.2 presents the reliability analysis. 10

Confidence regions

When uncertainties are taken into account, the graphic displayed in Fig. 3 (related to the deterministic analysis) becomes Figs. 5 and6 for the two values 0.05 and 0.10 of δ, the parameter that controls the level of uncertainty. Note that only one half of the eigenvalues are plotted in order to get a better visualization. Figures 5 and6 show the mean and the 95% confidence region of the corresponding real and imaginary parts of the eigenvalues. When δ increases, the confidence region also increases, as expected. Also expected is the fact that when u = 0 there is no uncertainty in the random eigenvalues, since for u = 0 the fluid model does not play any role. However, one thing difficult to predict is how the robustness of the system varies with speed u. For example, Fig. 6 shows that, for u between 7 and 8, the confidence region of Re(Λ nd ) and Im(Λ nd ) is narrower and close to u = 12 it is broader.

An interesting aspect of the stochastic response for small values of u. If we take a closer look on the probabilistic envelopes of the real part of the eigenvalues depicted in Fig. 7. It can be noted that for different values of the reference frequency, the real part of the eigenvalues are different from zero, and the influence region depends on this frequency. For lower frequencies, this region is smaller. When ω r is higher, the entries of the nominal value of the stochastic operator is higher; see Eq. ( 15). Since the same δ is used, the dispersion of the matrix norm will be different depending on the value of ω r . Consequently, it has a different impact on the stochastic response. This means that the reference frequency is an additional hyperparameter in the stochastic modeling, which can either be identified if experimental data are available or be used as a sensitivity parameter in a robust analysis with respect to model uncertainties if not. Either way, the conclusion is that this result is valuable to understand that modeling errors might influence the real part of the eigenvalues in low speeds, which should be a concern when analyzing stability.

In the frequency domain the response shown in Fig. 4 becomes random when uncertainties are taken into account. Figure 8 shows the response for u = 2.5 and δ = 0.1, where the gray region represents the 95% confidence region. To better visualize the uncertainties, Fig. 9 shows how the coefficient of variation (standard deviation divided by the mean) changes with the frequency. It is noted that the level of uncertainty is higher close to the natural frequencies. This should be a concern for the designer of this system, depending on the excitation frequencies.

These results are valuable to investigate the propagation of uncertainty throughout the computational model. It is not evident, a priori, how uncertainties related to the fluid-structure interaction will impact the eigenvalues of the system. With parameter β fixed it is possible to analyze the uncertainty of the response for different parameters of δ and u. When the confidence limit is small, the model is robust to global uncertainties related the fluid-structure interaction. On the other hand, if the confidence limit is big, the model response is not reliable.

Reliability analysis

When the lines of the graphics cross each other (Fig. 6, for example), it is difficult to follow the different confidence regions. Instead of analyzing these crossings, let us evaluate the probability of the occurrence of instability. This probability of occurrence of instability is the failure probability of our system.

Two probabilities are defined: (1) one related to the occurrence of divergence instability mode and (2) other related to the flutter instability mode. Figure 10(a) shows that until about u = 2.5 the probability of divergence is close to zero, but for u > 4 it becomes almost 100%. At the same time, Figure 10(b) shows that until about u = 5 the probability of flutter is close to zero, but for u > 7.5 it becomes almost 100%.

Figure 11 shows the estimation of the probability density function of Re(Λ 2nd ) for two different values of u. In both cases there is an important probability that Re(Λ 2nd ) be in the neighborhood of zero. However, for higher flow speeds, the probability that this random variable be in a region greater than zero increases. An interesting fact to note is that random variable Re(Λ 2nd ) follows a mixed probability density function, with a concentrated probability in the neighborhood of zero (stable condition, if only flutter instability is considered) and a continuous probability density function above zero (unstable condition related to flutter). The concentrated probability mass near zero approaches the delta Dirac function as the number of Monte Carlo simulations approaches infinity. Now we will detail the analysis of the probability of being in a flutter condition for different levels of uncertainties. Figure 12 shows the probability of flutter as a function of the flow speed u for different levels of uncertainty (delta values). The solid arrows show that if the deterministic (δ = 0) limit speed is considered, the probability of occurrence of instability (or failure probability) would be about 50%, which means that uncertainties should be taken into account.

On the other hand, the dashed-dot arrow shows that if we define a failure probability, say 1%, depending on the level of uncertainty, then there will be a limit flow speed for each value of δ: u ≃ 5.85 for δ = 0.1, u ≃ 6.06 for δ = 0.05 and u ≃ 6.24 for δ = 0.01. As expected, if δ increases, the limit u should decrease because the scenario is getting worse.

Finally, we construct a control surface that might help a designer to make a decision based on the values of the fluid speed u and on the level δ of uncertainty. This surface, Fig. 13, shows the probability of being in a flutter region for different values of u and δ: u ∈ [5.6, 7.0] and δ ∈ [0.0, 0.1]. Note that fixing the value of δ,a su increases, the probability of flutter always increases. On the other hand, it is interesting to note that fixing the value of u,a sδ increases the probability of flutter might increase or decrease. Either way, depending on the level of uncertainty and on the value of the fluid speed, one knows the probability of being unstable.

Although the results of the previous section are valuable to investigate the propagation of uncertainty, we would like to use the computation model for design. In this sense, we are interest in avoiding failures, which are related to instabilities in the present analysis. Instead of setting a safety factor, a more refined analysis (reliability analysis) is done with the stochastic computational model, where the probability of failure is computed. For a fixed β, the surface constructed in Fig. 13 is the tool needed by the designer to visualize the level of risk he is taking for each pair (δ,u).

Conclusions

The stability problem of a pipe conveying fluid, [START_REF] Paidoussis | Dynamic stability of pipes conveying fluid[END_REF], has been extended to the case for which uncertainties are taken into account in the computational model. A simple Euler-Bernoulli beam and a plug flow model are considered to model the dynamical system, which is discretized by means of the finite element methods and reduced using eigenmodes of the beam. Two analyzes are done: (1) construction of the confidence regions for the eigenvalues and for the stochastic response and (2) reliability analysis related to the probability of being in a divergence region or in a flutter region. Both analyzes have been done varying the dimensionless speed of the fluid and the level of uncertainty of the probabilistic model. An interesting result is that a mixed distribution has been obtained for the real part of one random eigenvalue.

A stochastic structure based on the nonparametric probabilistic approach of model uncertainties has been proposed for the operator related to the noninertial fluid forces. Therefore, the proposed stochastic model considers globally uncertainties in the noninertial fluid forces, but other uncertainties, such as the ones related to the mass properties, to the structure, or to the forces acting on the system are not taken into account. This approach is built upon the introduction of a stochastic structure modifying directly the operators that underline the physics of the model, in clear contrast with the parametric alternatives, in which uncertain parameters are modeled as random variables. This nonparametric formulation, which easily can be combined with parametric ones, is employed to capture structural model discrepancies, and, therefore, enhances the ability of making reliable predictions.

It is the first time such a stochastic structure is proposed, therefore, more analysis must be pursued. For example, the prediction model should be compared with sophisticated models or high fidelity models and experimental data.

List of captions: Figure 2. Venn diagram to illustrate the difference between parametric and nonparametric probabilistic approaches. Figure 7. For small values of u and δ = 0.1, 95% confidence limits of the real part of the eigenvalues.

Figure 8. Frequency response amplitude at η = 0.5 for u = 2.5 and δ = 0.1. Figure 9. Coefficient of variation of the frequency response amplitude at η = 0.5 for u = 2.5.

Figure 1 .

 1 Figure 1. Sketch of the system considered in the analysis (the arrow represents the internal fluid flow).Figure 2. Venn diagram to illustrate the difference between parametric and nonparametric probabilistic approaches. Figure 3. Stability charts. Real (a) and imaginary components (b) of λ nd . Figure 4. Frequency response amplitude at η = 0.5 for two different values of u. Figure 5. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.05. (a) real and (b) imaginary components. Figure 6. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.10. (a) real and (b) imaginary components.Figure7. For small values of u and δ = 0.1, 95% confidence limits of the real part of the eigenvalues.Figure8. Frequency response amplitude at η = 0.5 for u = 2.5 and δ = 0.1. Figure9. Coefficient of variation of the frequency response amplitude at η = 0.5 for u = 2.5. Figure10. For different values of u and δ = 0.1, (a) probability of occurrence of divergence instability and (b) probability of occurrence of flutter instability.Figure11. Approximation of the probability density function of Re(Λ 2nd ) for two different values of u.Figure12. For different values of u and δ, probability of occurrence of flutter instability. Figure13. Surface showing the probability of occurrence of flutter instability for different values of u and δ.

Figure 3 .

 3 Figure 1. Sketch of the system considered in the analysis (the arrow represents the internal fluid flow).Figure 2. Venn diagram to illustrate the difference between parametric and nonparametric probabilistic approaches. Figure 3. Stability charts. Real (a) and imaginary components (b) of λ nd . Figure 4. Frequency response amplitude at η = 0.5 for two different values of u. Figure 5. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.05. (a) real and (b) imaginary components. Figure 6. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.10. (a) real and (b) imaginary components.Figure7. For small values of u and δ = 0.1, 95% confidence limits of the real part of the eigenvalues.Figure8. Frequency response amplitude at η = 0.5 for u = 2.5 and δ = 0.1. Figure9. Coefficient of variation of the frequency response amplitude at η = 0.5 for u = 2.5. Figure10. For different values of u and δ = 0.1, (a) probability of occurrence of divergence instability and (b) probability of occurrence of flutter instability.Figure11. Approximation of the probability density function of Re(Λ 2nd ) for two different values of u.Figure12. For different values of u and δ, probability of occurrence of flutter instability. Figure13. Surface showing the probability of occurrence of flutter instability for different values of u and δ.

Figure 5 .

 5 Figure 1. Sketch of the system considered in the analysis (the arrow represents the internal fluid flow).Figure 2. Venn diagram to illustrate the difference between parametric and nonparametric probabilistic approaches. Figure 3. Stability charts. Real (a) and imaginary components (b) of λ nd . Figure 4. Frequency response amplitude at η = 0.5 for two different values of u. Figure 5. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.05. (a) real and (b) imaginary components. Figure 6. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.10. (a) real and (b) imaginary components.Figure7. For small values of u and δ = 0.1, 95% confidence limits of the real part of the eigenvalues.Figure8. Frequency response amplitude at η = 0.5 for u = 2.5 and δ = 0.1. Figure9. Coefficient of variation of the frequency response amplitude at η = 0.5 for u = 2.5. Figure10. For different values of u and δ = 0.1, (a) probability of occurrence of divergence instability and (b) probability of occurrence of flutter instability.Figure11. Approximation of the probability density function of Re(Λ 2nd ) for two different values of u.Figure12. For different values of u and δ, probability of occurrence of flutter instability. Figure13. Surface showing the probability of occurrence of flutter instability for different values of u and δ.

Figure 6 .

 6 Figure 1. Sketch of the system considered in the analysis (the arrow represents the internal fluid flow).Figure 2. Venn diagram to illustrate the difference between parametric and nonparametric probabilistic approaches. Figure 3. Stability charts. Real (a) and imaginary components (b) of λ nd . Figure 4. Frequency response amplitude at η = 0.5 for two different values of u. Figure 5. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.05. (a) real and (b) imaginary components. Figure 6. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.10. (a) real and (b) imaginary components.Figure7. For small values of u and δ = 0.1, 95% confidence limits of the real part of the eigenvalues.Figure8. Frequency response amplitude at η = 0.5 for u = 2.5 and δ = 0.1. Figure9. Coefficient of variation of the frequency response amplitude at η = 0.5 for u = 2.5. Figure10. For different values of u and δ = 0.1, (a) probability of occurrence of divergence instability and (b) probability of occurrence of flutter instability.Figure11. Approximation of the probability density function of Re(Λ 2nd ) for two different values of u.Figure12. For different values of u and δ, probability of occurrence of flutter instability. Figure13. Surface showing the probability of occurrence of flutter instability for different values of u and δ.

Figure 10 .

 10 Figure 1. Sketch of the system considered in the analysis (the arrow represents the internal fluid flow).Figure 2. Venn diagram to illustrate the difference between parametric and nonparametric probabilistic approaches. Figure 3. Stability charts. Real (a) and imaginary components (b) of λ nd . Figure 4. Frequency response amplitude at η = 0.5 for two different values of u. Figure 5. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.05. (a) real and (b) imaginary components. Figure 6. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.10. (a) real and (b) imaginary components.Figure7. For small values of u and δ = 0.1, 95% confidence limits of the real part of the eigenvalues.Figure8. Frequency response amplitude at η = 0.5 for u = 2.5 and δ = 0.1. Figure9. Coefficient of variation of the frequency response amplitude at η = 0.5 for u = 2.5. Figure10. For different values of u and δ = 0.1, (a) probability of occurrence of divergence instability and (b) probability of occurrence of flutter instability.Figure11. Approximation of the probability density function of Re(Λ 2nd ) for two different values of u.Figure12. For different values of u and δ, probability of occurrence of flutter instability. Figure13. Surface showing the probability of occurrence of flutter instability for different values of u and δ.

Figure 12 .

 12 Figure 1. Sketch of the system considered in the analysis (the arrow represents the internal fluid flow).Figure 2. Venn diagram to illustrate the difference between parametric and nonparametric probabilistic approaches. Figure 3. Stability charts. Real (a) and imaginary components (b) of λ nd . Figure 4. Frequency response amplitude at η = 0.5 for two different values of u. Figure 5. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.05. (a) real and (b) imaginary components. Figure 6. Mean and 95% confidence interval of one half of Λ nd 's for δ = 0.10. (a) real and (b) imaginary components.Figure7. For small values of u and δ = 0.1, 95% confidence limits of the real part of the eigenvalues.Figure8. Frequency response amplitude at η = 0.5 for u = 2.5 and δ = 0.1. Figure9. Coefficient of variation of the frequency response amplitude at η = 0.5 for u = 2.5. Figure10. For different values of u and δ = 0.1, (a) probability of occurrence of divergence instability and (b) probability of occurrence of flutter instability.Figure11. Approximation of the probability density function of Re(Λ 2nd ) for two different values of u.Figure12. For different values of u and δ, probability of occurrence of flutter instability. Figure13. Surface showing the probability of occurrence of flutter instability for different values of u and δ.

Figure 1 :

 1 Figure 1: Sketch of the system considered in the analysis (the arrow represents the internal fluid flow).

Figure 2 :Figure 3 :

 23 Figure 2: Venn diagram to illustrate the difference between parametric and nonparametric probabilistic approaches.

Figure 4 :

 4 Figure 4: Frequency response amplitude at η = 0.5 for two different values of u.

Figure 7 :Figure 8 :

 78 Figure 7: For small values of u and δ = 0.1, 95% confidence limits of the real part of the eigenvalues.

Figure 9 :Figure 10 :

 910 Figure 9: Coefficient of variation of the frequency response amplitude at η = 0.5 for u = 2.5 and δ = 0.1.

Figure 11 :

 11 Figure 11: Approximation of the probability density function of Re(Λ 2nd ) for two different values of u.

Figure 12 :Figure 13 :

 1213 Figure 12: For different values of u and δ, probability of occurrence of flutter instability .

Acknowledgements

The authors are grateful to the financial support from the Brazilian agencies: CAPES, CNPq, and FAPERJ. And specifically the project CAPES-COFECUB.


Sébastien Cargou

Marc Fouet

Armelle Montrose

David Bourrier

et al.. AN INTE- GRATED MAGNETIC PLANAR Karine Reybier

Anne Marie Gué

AN INTEGRATED MAGNETIC PLANAR ACTUATOR REDEFINING MULTILEVEL (3D) MICROFLUIDIC STRATEGIES

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

AN INTEGRATED MAGNETIC PLANAR ACTUATOR REDEFINING MULTILEVEL (3D) MICROFLUIDIC STRATEGIES

S. CARGOU 13 , M. FOUET 13* , A. MONTROSE 23 , D. BOURRIER 13 , K. REYBIER 23 , A.-M. GUÉ 13 1 LAAS-CNRS, Toulouse, FRANCE 2 PHARMADEV, Toulouse, FRANCE 3 Université de Toulouse, Toulouse, FRANCE *email : mfouet@laas.fr Super-paramagnetic micro particles (SMPs) are broadly used in medical and biological applications (from cellular to genomic scales), but the handling tools for those particles remain mostly non-integrated. Since a decade, a few solutions for integrated actuation have been proposed and paved the way to the micro-scale control of SMPs. So far, spiral planar coils showed the best compromises between the available magnetic force and the actuator footprint. Albeit the technology to integrate this type of coils with microfluidic structures exists, it is hard to find examples of dual architectures using those two aspects for concrete applications.

As a matter of fact, the fluidic structures of such an actuator should be designed so as to take advantage of the high impact of the magnetic force component developed on the orthogonal axis of a spiral planar micro-coil, which is impossible given the lack of simple fabrication techniques for 3D microdluidics.

We therefore developed a SU-8 dry film lamination technique, which can be used with standard micro technology processes. It allows us to combine a good resolution (the height and width of the coil turns are 5µm) and 3D structures, can be carried out on silicon, glass or flexible substrates, and is compatible with other commercially available low cost dry films for the microfluidic architecture. Our approach was then to place planar coils under two superposed channels, the top one containing a sample with magnetically targeted cells, and the bottom one a buffer to recover the magnetic beads an attached cells.

The devices were then tested with a setup including a current source (100mA injected in each coil, the coils being wired in three series of three), a Fluigent pressure controller, and a bright field microscope. In a first step (A), the magnetic particles (5µm diameter from Invitrogen) are vertically focused using the first six coils. Once the beads have reached the separation level, they are attracted towards the bottom channel with the last three actuators.

Using that system, we achieved to separate a solution of SMPs with an 80% efficiency at a 2.5 µL/min flow rate (4.6 mm/s). It is also possible to perform trapping (all the beads stay on the coils) at lower flow rate: trapping is observed simultaneously with separation (with a 100% yield) at 1µL/min (1.8 mm/s). We showed also that trapping of monocytes was possible on this coil design.

Currently we are working on the system characterization with a biological application. This type of systems is particularly suited to integrated detection techniques, and a next important step would be the conception of a fully coupled magnetic actuation and electrical detection microfluidic chip.

Figure 2 :Figure1:

 2 Figure 2: a) Scheme of the integrated magneto-fluidic actuator. SMPs are separated from the top channel to be released in the bottom channel, due to the magnetic attraction of nine microcoils. b) Top view photography of the microfluidic chip.