V Ammovilli
email: vammovilli@unistra.fr

M Bilasse

I Charpentier

Continuous nonlinear eigenvalue solver with applications to the design of electro/magnetorheological sandwich structures

Keywords: Smart sandwich structures, nonlinear free vibration, homotopy method, Taylor series, automatic differentiation

Smart sandwich structures comprising an electro-or a magnetorheological material have the potential to attenuate vibration over a wide range of frequencies. The analysis of their vibration behaviour with respect to the continuous variation of the field intensity is thus a major challenge for research and industry to maximize damping treatments. The numerical higher order homotopy method we propose models the effects of a continuous variation of the field intensity on resonant frequencies and loss factors by means of Taylor expansions. Comparisons between our continuous approach and the classical incremental method are proposed for state of the art sandwich beams and plate structures comprising ER/MR fluids to highlight the benefits of our continuous methods in terms of maximal damping determination.

Introduction

Sandwich structures are defined as multilayered composite structures, optimized for desired/anticipated lifetime loading conditions occurring in engineering applications such as automotive, railway, submarine, and aerospace. The analysis of their vibration behaviour is a major challenge for industry and research from experimental and numerical points of view.

The benefits of sandwich structures involving viscoelastic (VE), electrorheological (ER) or magnetorheological (MR) damping treatments are demonstrated for vibration attenuation or suppression [START_REF] Nashif | Vibration damping[END_REF][START_REF] Yalcintas | [END_REF]. Sandwich structures with their fixed VE material parameters usually yield limited vibration attenuation within a narrow frequency band. Recent advances on VE sandwich structures are reviewed in [3].

The rheological properties (viscosity and shear modulus) of ER/MR fluids or elastomers are functions of the applied external field respectively. When such a material is embedded between two elastic layers [START_REF] Carlson | Electrorheological fluid composite structures[END_REF][START_REF] Weiss K D, Duclos T G, Chrzan | Magnetorheological fluid composite structures[END_REF], the stiffness and damping properties of the sandwich structure can be controlled by an external electric/magnetic field in a reversible manner. Sandwich structures with ER/MR material are known to offer attractive potentials for vibration attenuation over a wide range of frequencies and temperatures [START_REF] Narayana | [END_REF]7]. Moreover, additive manufacturing such as 3D printing allows for the design of complex core materials, honeycomb filled with a MR elastomer for instance [8], improving the mechanical stability, stiffness and energy absorption of the structure. Recent state-of-the-art reviews are discussed in [9][10][11].

As for VE models, classical ER/MR models combine elastic springs and viscous dashpots to describe the rheological properties. These can be solid-like models such as Kelvin-Voigt and Zener elements model, and fluid-models such as Maxwell fluid and three-parameter fluid model [12]. Rheological parameters are characterized from experiments carried out on ER/MR fluids or elastomers [9]. The shear modulus G * c can be written as a function of the driving frequency ω and the varying intensity p of the applied field over some interval I:

G * c (ω, p) = G c (ω, p) 1 + iη c (ω, p) , (1)
where the real part G c (ω, p) is the storage modulus and the ratio η c (ω, p) represents the loss factor of the core material. The imaginary part G c (ω, p)η c (ω, p) is the loss modulus. The rheological law (1) may be implemented in the finite element method (FEM) to predict the dynamic behaviour of a sandwich structure [12][13][14][15] or to study its sensitivity with respect to modelling parameters such as thickness ratio, material properties or field intensity values [16][17][18][19][20][21]. Note that most of the finite element analyses (FEA) assume a constant frequency shear modulus. When a frequency dependent core is embedded in a structure, the free vibration analysis based on FEM can be formulated as the generalised nonlinear eigenvalue problem:

K ω(p), p -ω 2 (p)M(p) u(p) = 0, ( 2
)
where K is the stiffness matrix, M is the mass matrix and (u(p),ω(p)) is the eigensolution. Structural damping properties (Ω(p), η(p)) are deduced following:

ω 2 (p) = Ω 2 (p) 1 + iη(p) , for p ∈ I, (3)
where Ω(p) is the resonant frequency and η(p) is the structural loss factor. For a prescribed parameter p, the nonlinear frequency dependent problem (2) can be solved as proposed in [22], for instance. This algorithm combines homotopy to a higher order continuation method and agrees with the VE rheological modelling. Intensity effects can be thus evaluated at given values p of the interval. The calculation of continuous approximations for (u(p), ω(p)) should take into account the continuous variation of the modelling parameter p ∈ I. A very few methods propose an actual continuous approach in their parametric analyses. Among them, He [23] presents some analytical perturbation and homotopy analytical techniques for the solution of nonlinear differential equations by scaling the modelling parameter p by the homotopy parameter. Likewise, Charpentier and Lampoh [24] combine FEM, homotopy and higher order continuation to analyse the sensitivity of eigensolutions with respect to a damping function scaled by the homotopy parameter. Building on previous works [22,[25][26][27], Akoussan et al [28,29] propose a higher order continuation method for the sensitivity analysis of the damping properties of frequency dependent VE composite multilayer plates with respect to geometrical modelling parameters such as the layer thickness ratio and the fibers' orientation, respectively.

The calculation of the vibration behaviour over a large frequency range and/or a large field intensity range is invaluable from a design point of view to determine maximum damping configurations. To the best of our knowledge, any method has been proposed to evaluate the damping capabilities over a continuous range of electric/magnetic field intensity. The homotopy-based numerical method we present builds on [22] to fill this gap.

The manuscript is organized as follow. Section 2 presents the nonlinear free vibration problem associated to the FEM. The two-step homotopy method proposed for the continuous vibration analysis is described with detail in section 3. Numerical results, carried out on ER/MR beams and plates published in the literature, are discussed in section 4. Section 5 provides concluding remarks and future perspectives.

Setting of the residual nonlinear problem

For the sake of generality [22], the nonlinear eigenvalue problem (2) is written in the residual form:

R u(p), λ(p), p = K(λ(p), p) -λ(p)M(p) u(p) = 0, (4)
where u(p) = u ∈ C N denotes the eigenmode, λ(p) = λ ∈ C is the eigenvalue such that λ = ω 2 , and p ∈ I = [p min , p max ] represents the modelling parameter. In the paper, the modelling parameter p to be varied continuously is the intensity of the electric or magnetic field.

The residual problem (4) can be written in a simplified manner as:

R(u, λ, p) = K(λ, p) -λM u = 0, (5)
since the mass matrix of the ER/MR sandwich structure does not depend on the applied field. The stiffness matrix K(λ, p) is decomposed into:

K(λ, p) = K * 0 + G * (λ, p)K v , (6)
where K * 0 = K(0, 0) is the finite element delayed elasticity stiffness matrix computed assuming λ = 0 and p = 0, K v is a real matrix and G * (λ, p) is a complex nonlinear scalar function describing the dependence of the stiffness matrix in both the frequency and the field intensity. Introducing (6) into (5) yields the nonlinear general eigenvalue problem:

R(u, λ, p) = K * 0 + G * (λ, p)K v -λM u = 0. (7)
This nonlinear problem cannot be solved in a straightforward manner, even for a prescribed value of the field intensity p. The reader is referred to [30] for a review of nonlinear eigenvalue problems.

Generalised homotopy

Former studies [22,31] propose a continuous method for the solution of nonlinear eigenvalue problems. The main idea is to split the residual problem (7) into an eigenvalue problem S(u, λ) = 0 with known solutions and some additional nonlinear contribution T (u, λ, p), that is:

R(u, λ, p) = S(u, λ) + T (u, λ, p) = 0, (8)
where the functions S and T are specific to the problem under study (see table 1) and case studies presented in section 4.

As described in subsection 3.1, the solutions of S(u, λ) = 0 can be then deformed to compute the solutions of R(u, λ) = 0 through the so-called homotopy method.

The generalized higher order homotopy method we propose moreover allows for the continuous solution of the nonlinear eigenvalue problem while conserving its mechanical/physical meaning. In others words, the continuous variations of some material or geometrical or intensity parameter p onto the eigenvalues (u(p), λ(p)), resonant frequencies Ω(p) or loss factors η(p) can be studied in an actual continuous manner.

The use and benefits of this new generalized homotopy-based eigensolver are exemplified in subsection 3.2 and section 4.

General nonlinear eigensolver

The solution of ( 8) is implemented through a homotopy method mapping S to R by computing continuous branches of solutions for: R u(α), λ(α), p(α) = S(u, λ) + αT (u, λ, p) = 0, (9) for the homotopy parameter α ∈ R varying in the range [α min , α max ].

The unknowns are developed into Taylor series in the vicinity of α to be able to evaluate continuous approximation branches:

λ(α + δ α , p) ≈ N n=0 (δ α ) n n! ∂ n λ ∂α n (α, p) = N n=0 (δ α ) n λ n , (10)
u(α + δ α , p) ≈ N n=0 (δ α ) n n! ∂ n u ∂α n (α, p) = N n=0 (δ α ) n u n , (11)
where δ α is a small perturbation of α and N is the truncation order. Unknowns u n and λ n are the Taylor coefficients of u and λ, respectively. As noted in [23,24], the linear variation of the modelling parameter p can be accounted by scaling it with respect to the homotopy parameter α ∈ R. However, see section 4, most of the ER/MR constitutive laws cannot be reformulated to exhibit such a linear dependence.

Any nonlinear bijective function mapping p to α, for instance F (p(α)) = α and F -1 (α) = p, can be accounted through a higher order differentiation of F (p) with respect to the homotopy parameter α. This is implemented by means of the chain rule and the Taylor expansion of p(α) with respect to α:

p(α + δ α ) ≈ N n=0 (δ α ) n n! ∂ n p ∂α n (α) = N n=0 (δ α ) n p n , (12)
where p n are the Taylor coefficients of p. From a practical point of view, the Taylor series (12) implements the variation of ER/MR properties parametrized by p. The higher order approximations ( 10)-( 12) then allow for the calculation of continuous branches of solutions with respect to the variation of p in the range [p min , p max ].

Introducing these series into (9) yields:

R u(α + δ α ), λ(α + δ α ), p(α + δ α ) ≈ N n=0 (δ α ) n R n = 0, ( 13
)
where R n are the Taylor coefficients of R. The analyticity assumption implies:

R n = S n + αT n = 0 for n = 1, .., N, (14)
where the Taylor coefficients S n and T n of functions S and T satisfy:

S u(α + δ α ), λ(α + δ α ) ≈ N n=0 (δ α ) n S n , (15)
and

T u(α + δ α ), λ(α + δ α ), p(α + δ α ) ≈ N n=0 (δ α ) n T n . (16)
Faà di Bruno's generalised chain rule formula [22] allows to write a sequence of N + 1 equations for the unknown Taylor coefficients u n and λ n .

At order n = 0, equation ( 14) is restricted to:

J u 0 = S 1 | u1=Id , λ1=0 + α (j) T 1 | u1=Id , λ1=0, p1=0 u 0 = 0, (17)
where J , {S 1 | u1=Id , λ1=0 } and {T 1 | u1=Id , λ1=0 } are Jacobian matrices, differentiated with respect to u. One deduces that t u 0 J = 0 since the finite element matrices contained in J are symmetric. At order n, equation ( 14) is split into some first order terms multiplied by the unknown Taylor coefficients u n and λ n and some complementary higher order terms. One deduces:

S 1 | u1=Id , λ1=0 + α (j) T 1 | u1=Id , λ1=0, p1=0 u n + (18) + S 1 | u1=0 , λ1=1 + α (j) T 1 | u1=0 , λ1=1, p1=0 λ n + (19) + S n | un=0 , λn=0 -α (j) T n | un=0 , λn=0 -T n-1 = 0. ( 20
)
For the sake of clarity, ( 20) is rewritten as:

J u n + Dλ n + H n = 0, (21)
where the gradient vector D, differentiated with respect to u and λ, satisfies:

D = S 1 | u1=0 , λ1=1 + α T 1 | u1=0 , λ1=1 , p1=0 . (22)
The nonlinear higher order contribution at order n is:

H n = S n | un=0 , λn=0 + α T n | un=0 , λn=0 + T n-1 . (23)
Multiplying ( 21) by t u 0 yields the following equation for the unknown λ n :

λ n = - t u 0 H n t u 0 D . ( 24
)
Equation ( 21) cannot be solved in a straightforward manner since the Jacobian J is singular. The orthogonality condition ( 25) is added to close the system:

t u 0 u(α) -u 0 = 0. ( 25
)
It implies that the eigenmode u 0 is orthogonal to the Taylor coefficients u n of u for n = 1, .., N .

A Lagrange multiplier κ is introduced [31] to compute u n following:

  J u 0 t u 0 0      u n κ    =    -Dλ n -H n 0    . (26)
The coefficients of the series ( 10) and ( 11) are obtained by computing the higher order term H n , by evaluating (24) and by solving the linear system of (26), in an alternate manner.

The convergence radius α * of the computed series is estimated as:

α * = ε u 1 u N 1 N -1 , ( 27
)
where ε is a small parameter (usually ε = 10 -6 ). A branch of solutions is then built from Taylor expansions of u and λ by evaluating the series up to α * .

For each mode, one may deduce Taylor expansions for the resonant frequency Ω 2 (p) and the loss factor η(p) thanks to the chain rule since they satisfy (3). Two important benefits can be highlighted. Firstly, continuous branches of solutions (p, Ω 2 (p)) and (p, η(p)) can be plotted as presented in section 4 for different structures and constitutive laws. Secondly, Taylor expansions can be used to determine minimum or maximum values in a range in a very efficient and accurate manner.

The computation of an additional solution branch is required as far as α max is not reached. The new branch is computed from the solution of the previous branch by setting (u 0 , λ 0 ) = u(α * ), λ(α * ) .

Note that all the derivative computation are carried out by means of automatic differentiation. [32]. This provides generality to the nonlinear continuous eigensolver.

Linear VE Problem

[K * 0 -λM] u = 0 2. Nonlinear VE Problem [K * 0 + αG * 1 (λ, p min )K v -λM] u = 0 3. Nonlinear ER/MR Problem [K * 0 + αG * 2 (λ, p)K v -λM] u = 0
S(u, λ) K * 0 -λM u K * 0 + G * 1 (λ, p min )Kv -λM u T (u, λ, p) G * 1 (λ, p min )Kvu G * 2 (λ, p)Kvu α α ∈ R F (p) p p min F -1 (α) = N n=0 (δα) n pn α min 0 F (p min ) αmax 1 1 G * (λ, p) G * 1 (λ, p min ) G * 2 (λ, p)

Specialisation for ER/MR law

The calculation of the vibration behaviour of an ER/MR sandwich structure over a continuous range [p min , p max ] of electric/magnetic field intensity can be carried out as follows. The general homotopy process described in subsection 3.1 is now specialised to agree with the algorithm described in figure 1. Therein, the homotopy solver is run twice from input data (parameters, parameter ranges, rheological laws and functions) reported in table 1.

The first homotopy is the one described in [22]. It deforms continuously the solution of the linear VE problem:

[K * 0 -λM] u = 0, (28)
related to α = 0 and computable using a classical eigensolver, to the solution for the nonlinear VE problem:

[K * 0 + αG * 1 (λ, p)K v -λM] u = 0, for α ∈ [0, 1], (29)
where the parameter p is temporarily set equal to p min , as written in table 1. However, in the vibration analysis of a VE sandwich structure [22], one is interested in this last solution only.

The second homotopy is devoted to the evaluation of the continuous variation of the eigensolution with respect to the modelling parameter p. Its input data are chosen to map the field intensity parameter p ∈ [p min , p max ] to the homotopy parameter α through the bijective scaling function α = F (p). The homotopy parameter thus varies in [F (p min ), 1].

Any actual continuous rheological law G * 1 (λ, p) can be accounted. It suffices to split it into some real-valued bijective scaling function F (p) and a complex-valued virtual law G * 2 (λ, p) satisfying:

G * 1 (λ, p) = F (p) G * 2 (λ, p). ( 30
)
The scaling function F (p) ∈ R can be identified to the homotopy parameter α ∈ R as far as it satisfies F (p max ) = α max = 1.

For efficiency reasons, it is convenient to choose F so that its reciprocal function F -1 is of simple computation. Where possible, one may also choose G * 2 independent from p. Typically, depending on the law, one may choose either F (p) = p/p max or F (p) as a polynomial fraction of p. Note that the choice of the function F may impact the number of iterations and the distance between two consecutive iterations.

When the constitutive law is independent from the frequency, it is worth mentioning that the first homotopy can be replaced by a complex eigenvalue solver applied to the linear VE problem.

Numerical results

The vibration analysis is carried out on three ER or MR sandwich structures described in the literature to assess the benefits of our continuous homotopy approach over the classical incremental approaches.

From a numerical point of view, the eigenvalue problems are discretized using the finite element method. The number of elements is specified for each case study in the tables describing the geometrical and material properties. The truncature order of the Taylor expansions is N = 20. With regards to formula (27), the convergence radius a max and the parameter ε are both set equal to 10 -6 . These data allow for reproducibility.

4.1.

Rectangular "Aluminium/ER fluid/Aluminium" sandwich plate with a frequency independent constitutive law The rectangular sandwich plate [33] presented in figure 2 comprises a lower host plate (layer 3) and an upper constraining plate (layer 1) made of aluminium alloy, and a core (layer 2) composed of an outer adhesive joint embedding an ER fluid. This study is

L x L y h 3 h 2 h 1 Figure 2.
Rectangular "Aluminium/ER fluid/Aluminium" sandwich plate described in [33].

based on the four-node plate element used in [34]. In the FEM, the core properties are assumed to be homogeneous. This plate is simply supported on all its edges (SSSS).

The geometrical and material properties of the structure are reported in table 2. The shear modulus of the core material G * c (E), as measured in [START_REF] Don | Theses and Dissertations[END_REF], depends on the electric field intensity E ranging from 0 kV.mm -1 to 3.5 kV.mm -1 :

G * c (E) = 15000E 2 + i6900, (31)
G * c (E) = G * 0 + G * 1 (E). (32)
As noted in subsection 3.2, the first homotopy is useless since the constitutive law (31) does not depend on frequency.

For each eigenmode of the VE problem, the continuous variation of the nonlinear eigenvalue with respect to the field intensity is operated by mean of the "second" homotopy.

The law (31) is here split into a constant part G * 0 = G * c (0) = i 6900 and a nonlinear part G * 1 (E) = 15000E 2 . Following the discussion proposed in subsection 3.2, the nonlinear part of (32) is written as:

G * 1 (E) = F (E) G * 2 (E) = α G * 2 (E). ( 33
)
The scaling function F (E) = α can be chosen as:

F (E) = E 2 E 2 max = α, (34)
then G * 2 (E) is equal to 15000E 2 max . Even other choices are possible for the function F (E), ( 34) is the simplest one that describes the range of interest [0 kV.mm -1 , 3.5 kV.mm -1 ] since it yields the independence of G * 2 with respect to both the frequency and the electric field.

The nonlinear continuous branches of the first four natural frequencies and modal loss factors of the ER rectangular sandwich plate are plotted in figure 3 for a comparison with the discrete natural frequencies and modal loss factors computed by the incremental approach used in [33]. The agreement at the discrete points obtained with the incremental method is very good. Minor discrepancies (less than 6%) can be observed in table 3 for both resonance frequencies and loss factors. These result from different finite element discretization. Moreover, in the kinematic assumptions, the normal stress is accounted in the core layer in our FEM while it is neglected in the FEM implemented in [33].

Table 2. Geometrical and material properties of the rectangular "Aluminium/ER fluid/Aluminium" plate described in [33]. Note that points where series were computed are marked with a × to indicate the number of iterations. One observes that 7 series computations are needed to plot the variation of the first natural frequency f with respect to the variation of the electric field intensity E.

Plate

Rubber sealed "Aluminium/MR fluid/Aluminium" sandwich beam with a frequency independent constitutive law

The rectangular sandwich beam [START_REF] Rajamohan | [END_REF] presented in figure 4 is studied under clamped-free (CF) boundary conditions. This study is based on the two-node beam element used in [37]. The external elastic face layers (layers 1 and 3) are made of aluminium alloy . Rubber sealed "Aluminium/MR fluid/Aluminium" sandwich beam described in [START_REF] Rajamohan | [END_REF].

and the core layer (layer 2) is made of a MR fluid embedded in a silicon rubber sealant. Their shear moduli are assumed to be independent from the frequency. The magnetic field B is varied from 0 mT to 50 mT. The geometrical and material properties of the MR sandwich beam are described in table 4.

As in [START_REF] Rajamohan | [END_REF], the core layer of the sandwich beam comprising the rubber seal and the MR fluid is modelled as a homogeneous material layer. The equivalent shear modulus is expressed as a combination of the moduli and widths of the two materials using the Table 4.

Geometrical and material properties of the rubber sealed "Aluminium/MR fluid/Aluminium" sandwich beam described in [START_REF] Rajamohan | [END_REF]. Table 5. Coefficients for the shear modulus formula modelling the viscoelastic behaviour of the MR fluid in the preyield regime (from [START_REF] Rajamohan | [END_REF]). so-called "rule of mixture":

Beam

G * c (B) = G r b r b + G * M R (B) 1 - b r b , (35)
where b r and b are the widths of the rubber and beam respectively, and G r and G * c are the shear modulus of the rubber and MR fluid respectively. In the preyield regime, the MR fluid exhibits a viscoelastic behaviour that can be modelled as [START_REF] Rajamohan | [END_REF]:

G * M R (B) = a 0 + a 1 B + a 2 B 2 + i b 0 + b 1 B + b 2 B 2 , (36)
where the coefficients a i and b i are provided in table 5.

The homotopy 1 (see figure 1) is applied by splitting [START_REF] Rajamohan | [END_REF] into:

G * c (B) = G * 0 + G * 1 (B), (37)
where

G * 0 = G * c (0) = G r b r b + (a 0 + ib 0 ) 1 - b r b , (38)
and

G * 1 (B) = B 1 - b r b a 1 + a 2 B + i(b 1 + b 2 B) . (39)
As discussed in subsection 3.2, the shear modulus (37) is split into:

G * 1 (B) = αG * 2 (B), ( 40
)
where α is chosen as a simple real-valued scaling function of B such that F (B max ) = 1: The complex-valued and field dependent function G * 2 (B) then satisfies:

α = F (B) = B B max . (41)
G * 2 (B) = B max 1 - b r b a 1 + a 2 B + i(b 1 + b 2 B) . (42)
Under these assumptions, the nonlinear eigenproblem ( 7) is quadratic [30] in the field intensity B. Figure 5 plots the nonlinear branches of resonant frequencies and modal loss factors computed in the full range of validity [0 mT, 50 mT] of the rheological law using continuous solid, dashed and dotted lines. Continuous branches of solutions presented in figure 5 reflect the nonlinearity of the constitutive law with respect to the magnetic field intensity. Note that a unique branch of solutions suffices to describe the specified range of intensity because the range of validity of the series is here larger than 50 mT Resonant frequencies and loss factors computed for 0, 10, 20 and 30 mT by means of the incremental approach [START_REF] Rajamohan | [END_REF] are plotted as well.

A decay in the loss factor of the first mode is observed for larger B. The incremental method does not exhibit this decay because too few points were reported in [START_REF] Rajamohan | [END_REF]. A dichotomy method applied to the Taylor series of the loss factor of the first mode then indicates that the maximal damping occurs near to B 26.17 mT.

Numerical results are presented in table 6. A very good agreement with an error less than 2% is obtained between the resonant frequencies. A good agreement and similar trends are observed for the loss factors. One observes that the discrepancies on the loss factors are of the same order for the first three modes. These are mainly due to differences in the kinematic modelling. For instance, our model does not include the shearing effect in the face layers. The annular sandwich plate [9] presented in figure 6 comprises two external elastic layers made of aluminium alloy (layers 1 and 3) and a MR fluid core enclosed in a silicon rubber sealant (layer 2). The geometrical and material properties of the structure are reported in table 7. The annular plate is free at its inner edge and simply supported at its outer edge (FS). This study is based on the four-node plate element used in [34] with a radial mesh.

Table 7.

Geometrical and material properties of the rubber sealed "Aluminium/MRF 132DG/Aluminium" annular sandwich plate described in [9].

Plate

Outer radius Re = 250 mm Inner radius R i = 50 mm Finite element mesh 16×40 elements, radial mesh

Elastic faces

Young modulus E f = 6. In [9], the rubber and the MR fluid are homogenized as a single material. In the present FEM, we decide to implement them as two different materials.

On the one hand, the shear modulus of the MR fluid is modelled by a frequency and magnetic field dependent constitutive law [9]:

G * M R (f, B) = G (f, B) + iG (f, B), (43)
with

G (f, B) = a 0 + a 1 B + a 2 B 2 1 -e -a3f , (44)
and

G (f, B) = b 0 + b 1 B + b 2 B 2 1 -e -b3f , (45)
where the frequency f is in Hz and the magnetic field intensity B varies in [0 mT, 90 mT]. The coefficients a i and b i are provided in table 8. On the other hand, the rubber seal is modelled by a viscoelastic material, the shear modulus of which is described by [9]:

G r (f ) = 1340000 + i 152511 + 68.31f + 0.475f 2 , ( 46
)
where f is the frequency in Hz. Note that, the variations of the shear modulus regarding the frequency are negligible in the range of frequencies from 0 Hz to 100 Hz, thus G r can be approximated by the constant law:

G r = 1340000 + i152511. (47)
To apply the proposed homotopy method, the frequency and magnetic field dependent constitutive law (43) is split into:

G * M R (f, B) = G * 0 + G * 1 (f, B), (48)
where

G * 0 = G * M R (0, 0) = 0 and G * 1 (f, B) = G * M R (f, B).
Referring to section 3, (48), we propose the scaling:

G * 1 (f, B) = αG * 2 (f, B), (49)
where the parameter α of the second homotopy is chosen as:

α = F (B) = a 0 + a 1 B + a 2 B 2 a 0 + a 1 B max + a 2 B 2 max . (50)
The second homotopy method is then run considering:

G * 2 (f, B) = κ 1 -e -a3f + ig(B) 1 -e -b3f , (51)
where

κ = a 0 + a 1 B max + a 2 B 2 max , (52)
and

g(B) = b 0 + b 1 B + b 2 B 2 a 0 + a 1 B + a 2 B 2 . ( 53
)
Under these assumptions, ( 7) is a rational eigenvalue problem [30,38]. Computed frequencies and loss factors are reported in table 9. Results obtained by the double homotopy method (continuous variation of B) are in a very good agreement with those obtained by the incremental method (prescribed values for B, [9]). The differences in the loss factors result from the differences in the finite element modelling of the heterogeneous core properties.

Continuous branches of solutions presented in figure 7 reflect the nonlinearity of the constitutive law with respect to the magnetic field intensity. In particular, a decay of the loss factor of the first mode is observed for larger B. One notices that the incremental method does not exhibit this decay because too few points were reported in [9]. A dichotomy method applied to the Taylor series of the loss factor of the first mode then indicates that the maximal damping occurs near to B 73.43 mT.

Conclusion

The vibration analysis of smart rheological sandwich structures is carried out through a combination of numerical methods [22]. The finite element method is applied to the free vibration analysis in order to write the complex nonlinear eigenvalue problem in a matrix form. The homotopy method allows for the solution of such a complex eigenvalue problem. In addition, the continuous nature of any homotopy can be accounted through Taylor series approximations and automatic differentiation that provides generality to the proposed nonlinear continuous eigensolver. These methods allow for the implementation of a numerical eigenvalue solver that notably allows for to the analysis of the continuous variation of the resonant frequencies and the loss factors with respect to the continuous variation of one of the modelling parameters.

The proposed method is exemplified on three different ER and MR sandwich structures, considering the continuous variations of ER/MR laws proposed in the literature with respect to their field intensity parameter and/or the frequency. Calculation are carried out with respect to actual continuous variations of the intensity parameters. The continuous variations we present are in very good agreement with the related published discrete results. Where possible, we extended the computation of the resonant frequencies and loss factors over the whole range of validity of the constitutive law.

For a given mode, this continuous method allows for the exploration of the effects of the field intensity on the structural behaviour over a broad range. It notably highlights possible decays for larger field intensity variations, that were not exhibited in [9,[START_REF] Rajamohan | [END_REF] that use the incremental approach and a very few points to represent the loss factor variations. The maximal damping can then be determined by applying the dichotomy method to its Taylor series with little additional computational efforts.

E 1 = E 3 = 7 × 1 =Figure 3 .

 13713 Figure 3. Rectangular "Aluminium/ER fluid/Aluminium" sandwich plate. Comparison of the first four resonance frequencies and loss factors computed by means of the homotopy method (continuous solid/dashed or dotted lines) and the incremental methods (coloured markers).

Figure 4

 4 Figure 4. Rubber sealed "Aluminium/MR fluid/Aluminium" sandwich beam described in[START_REF] Rajamohan | [END_REF].

Gr = 1 . 34 × 10 6

 1346 Pa Width br = 1.5 mm Density ρr = 1233 kg.m -3

Figure 5 .

 5 Figure 5. Rubber sealed "Aluminium/MR fluid/Aluminium" sandwich beam. Comparison of the first three resonance frequencies and loss factors computed by means of the homotopy method (continuous solid/dashed or dotted lines) and the incremental methods (coloured markers).

5 ×Gr = 1 . 34 × 10 6

 51346 10 10 Pa Density ρ f = 2650 kg.m -3 Thickness h f Pa Width br = 2 mm Density ρ f = 1460 kg.m -3

Figure 7 .

 7 Figure 7. Rubber sealed "Aluminium/MRF 132DG/Aluminium" annular plate. Comparison of the first three resonance frequencies and loss factors computed by means of the homotopy method (continuous solid/dashed or dotted lines) and the incremental methods (coloured markers).

Table 1 .

 1 Figure 1. Structure of the general eigensolver for nonlinear continuous ER or MR eigenproblems. Input data for the homotopies carried out with respect to the frequency ω and the field intensity p, respectively.

		Homotopy 1
		(w.r.t. frequency ω)
		Homotopy 2
		(w.r.t. parameter p)
	Homotopy 1	Homotopy 2
	(w.

r.t. frequency ω) (w.r.t. parameter p)

Table 3 .

 3 Rectangular "Aluminium/ER fluid/Aluminium" sandwich plate. Comparison of the first four resonance frequencies, loss factors and relative errors.

	E (kV.mm -1 ) Homotopy		Incremental		Relative error
		method		method [36]
		f (rad.s -1 )	η	f (rad.s -1 )	η	f (%) η(%)
	0.5	159.00	0.0371 155.53	0.037 2.23	0.27
		354.98	0.0170 345.14	0.017 2.85	0.00
		446.39	0.0135 430.36	0.014 3.72	3.57
		641.72	0.0095 608.26	0.010 5.50	5.00
	1.5	170.38	0.0274 169.37	0.028 0.60	2.14
		367.48	0.0148 358.99	0.015 2.36	1.33
		459.02	0.0121 444.21	0.012 3.33	0.83
		654.63	0.0088 623.17	0.009 5.05	2.22
	3.5	203.90	0.0102 202.40	0.010 0.74	2.00
		415.04	0.0086 407.99	0.009 1.73	4.44
		509.81	0.0078 494.27	0.008 3.14	2.50
		709.39	0.0063 677.50	0.006 4.71	5.00
		L
		b		Aluminium
				h 1
				h 2
				h 3
			MR Fluid
	b r	Rubber

Table 6 .

 6 Rubber sealed "Aluminium/MR fluid/Aluminium" sandwich beam. Comparison of the first five resonance frequencies, loss factors and relative errors.

	B (mT) Homotopy		Incremental	Relative error
		method		method [36]
		f (Hz)	η	f (Hz)	η	f (%) η(%)
	0	41.67	0.0548 42.38	0.047 1.68	16.60
		227.42	0.0140 227.44 0.012 0.01	16.67
		620.44	0.0046 616.92 0.004 0.57	15.00
		1207.70 0.0022 1196.1 0.002 0.97	10.00
		1990.22 0.0013 1962.9 0.001 1.39	30.00
	10	44.24	0.0634 44.78	0.056 1.21	13.21
		231.14	0.0184 230.95 0.017 0.08	8.24
		623.85	0.0063 620.15 0.006 0.60	5.00
		1210.93 0.0031 1199.1 0.003 0.99	3.33
		1993.39 0.0019 1965.9 0.002 1.40	5.00
	20	46.22	0.0669 46.65	0.059 0.92	13.39
		234.22	0.0213 233.87 0.019 0.15	12.11
		626.74	0.0076 622.88 0.007 0.62	8.57
		1213.73 0.0038 1201.8 0.003 0.99	26.67
		1996.09 0.0022 1968.5 0.002 1.40	10.00
	30	47.73	0.0673 47.42	0.061 0.65	10.33
		236.70	0.0229 235.24 0.021 0.62	9.05
		629.12	0.0084 624.13 0.008 0.80	5.00
		1216.02 0.0042 1203.9 0.004 1.01	5.00
		1998.33 0.0025 1970.6 0.002 1.41	25.00
	h e
	h i
	h e
		R e	R i
				b r
	Figure 6. Rubber sealed "Aluminium/MRF 132DG/Aluminium" annular
	sandwich plate.
	4.3. Rubber sealed "Aluminium/MRF 132DG/Aluminium" annular sandwich plate
	with a frequency dependent constitutive law

Table 8 .

 8 Coefficients of the frequency and magnetic field dependent constitutive law (44)-(45) (from[9]).

	a 0	192161 Pa	b 0	45524.4 Pa
	a 1	30663.6 Pa.mT -1	b 1	812.4 Pa.mT -1
	a 2	243.625 Pa.mT -2	b 2	6.4412 Pa.mT -2
	a 3	0.00408 Hz -1	b 3	0.00074 Hz -1

Table 9 .

 9 Rubber sealed "Aluminium/MRF 132DG/Aluminium" annular plate. Comparison of the first three resonance frequencies, loss factors and relative errors.

	B ( mT) Homotopy	Incremental	Relative error
		method		method [9]
		f ( Hz)	η	f ( Hz)	η	f (%)	η(%)
	0	15.93	0.039 15.67	0.040 1.66	2.50
		43.14	0.021 41.96	0.017 2.81	23.53
		78.52	0.018 71.21	0.018 10.27 0.00
	30	18.74	0.124 18.67	0.123 0.37	0.81
		48.07	0.076 47.20	0.077 1.84	1.30
		87.91	0.073 85.79	0.071 2.47	2.82
	50	21.25	0.167 22.06	0.160 3.67	4.38
		52.38	0.109 52.09	0.105 0.56	3.81
		96.17	0.105 94.62	0.096 1.64	9.38
	90	27.13	0.182 28.87	0.174 6.03	4.60
		62.94	0.143 62.77	0.127 0.27	12.60
		116.31	0.135 115.29	0.119 0.88	13.45

The proposed homotopy based method is thus of great interest for the design and the vibration control of smart materials and structures. The same methodology can be applied to the analysis of sandwich structures with respect to temperature variations.


Wei Yeow

Pierre Choong

Maussion

Anne Laurent
email: laurent@lirmm.fr

Dominique Laurent

Summarizing Multidimensional Databases Using Fuzzy Rules

Keywords: Multidimensional Databases, Levelwise Algorithms, Fuzzy Rules, Data Visualization

In the context of multidimensional data, OLAP tools are appropriate for the navigation in the data, aiming at discovering pertinent and abstract knowledge. However, due to the size of the dataset, a systematic and exhaustive exploration is not feasible. Therefore, the problem is to design automatic tools to ease the navigation in the data and their visualization. In this paper, we present a novel approach allowing to build automatically blocks of similar values in a given data cube and to associate these blocks with rules. Our method is based on a levelwise algorithm (a la Apriori) and on the fuzzy set theory. The latter is considered here due to the fact that some of the blocks computed by our algorithm can overlap.

Introduction

Multidimensional databases have been studied and used for about 10 years [START_REF] Codd | Providing OLAP to user-analysts: An IT mandate[END_REF]. Relational databases have been intially proposed for On-Line Transactional Processing (OLTP) in order to facilitate and manage transactionoriented applications. However, relational databases are not suitable for the analysis of huge volumes of data, referred to as On-Line Analytical Processing (OLAP) in the literature. Multidimensional databases have been proposed to cope with this problem.

In this context, data are stored in multidimensional tables, the so-called data cubes, defined over several dimensions. The data within the cube is called measure. For instance, Figure 1 displays sale results in a data cube defined over two dimensions.

OLAP tools provide users with operations that ease manipulation of data in order to dicover relevant information.

However, data are often voluminous, and thus, an exhaustive exploration is impossible. Moreover, users of such systems are noncomputer scientists but rather decision makers or experts of the data. Although many reporting tools exist in commercial softwares, it is still very difficult to visualize multidimensional databases so that relevant information is automatically displayed. Cube dimensionality is normally greater than 4, which makes visualization difficult.

In this paper, an automatic processing is proposed to facilitate multidimensional data visualization. Our method identifies automatically sub-cubes of similar data. In the previous example, the sub-cube defined by products P 1, P 2 and city C1 is a block corresponding to the measure value 6. Since a block is a sub-cube, it can be associated with a rule. For instance, the rule associated with the block mentioned above is: If PRODUCT = P 1 or P 2 and CITY = C1 Then sales = 6.

In this case, all cells of the block contain the same measure value. However, it is not always the case. For instance, the block corresponding to products P 1, P 2, P 3 and cities C3 and C4 is mainly composed by the measure value 5. Therefore, one may consider it as being associated with 5. Similarly, a block associated with the measure value 2 can be considered for products P 3, P 4 and cities C4, C5, C6. These two blocks overlap since they share a cell for product P 3 and city C4. As illustrated by this example, blocks may overlap. Fuzzy logic is used in order to convey these overlappings, making it possible to represent information such as for product P 2 and in a lower manner for product P 3.

The goal of this work is to identify blocks and their overlappings (as shown in Figure 1) as quickly as possible and then to build the corresponing rules (fuzzy or not). The more the relevance of the blocks, the more the representation quality. In [START_REF] Choong | Computing Appropriate Representation for Multidimensional Data[END_REF], the authors study different manners to represent data cubes. In particular, the authors highlight the fact that some representations are more relevant than others since they allow to display the cube according to user-specified criteria. In this work, we consider as relevant the configurations where same measure values are grouped. Another organization is possible, as described in [START_REF] Choong | Computing Appropriate Representation for Multidimensional Data[END_REF], where the data are organized so that the measure is ordered in an increasing manner over all the dimensions. However, it is an NP-hard problem to automatically build such a representation of a cube. In this paper, the organization of data is not altered before the block identification, although it could be interesting.

¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¡ ¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ ¡¢ £ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ £ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ £ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ £ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ £ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ £ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¡£ ¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¡¤ ¥ ¡¥ ¡¥ ¡¥ ¡¥ ¡¥ ¥ ¡¥ ¡¥ ¡¥ ¡¥ ¡¥ ¥ ¡¥ ¡¥ ¡¥ ¡¥ ¡¥ ¥ ¡¥ ¡¥ ¡¥ ¡¥ ¡¥ ¥ ¡¥ ¡¥ ¡¥ ¡¥ ¡¥ ¥ ¡¥ ¡¥ ¡¥ ¡¥ ¡¥ ¥ ¡¥ ¡¥ ¡¥ ¡¥ ¡¥ ¥ ¡¥ ¡¥ ¡¥ ¡¥ ¡¥ ¥ ¡¥ ¡¥ ¡¥ ¡¥ ¡¥ ¥ ¡¥ ¡¥ ¡¥ ¡¥ ¡¥ ¦ ¡¦ ¡¦ ¡¦ ¡¦ ¡¦ ¦ ¡¦ ¡¦ ¡¦ ¡¦ ¡¦ ¦ ¡¦ ¡¦ ¡¦ ¡¦ ¡¦ ¦ ¡¦ ¡¦ ¡¦ ¡¦ ¡¦ ¦ ¡¦ ¡¦ ¡¦ ¡¦ ¡¦ ¦ ¡¦ ¡¦ ¡¦ ¡¦ ¡¦ ¦ ¡¦ ¡¦ ¡¦ ¡¦ ¡¦ ¦ ¡¦ ¡¦ ¡¦ ¡¦ ¡¦ ¦ ¡¦ ¡¦ ¡¦ ¡¦ ¡¦ ¦ ¡¦ ¡¦ ¡¦ ¡¦ ¡¦ § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § § ¡ § ¡ § ¡ § ¡¡¡¡ ¡¡¡ ¡¡¡ ¡¡¡ ¡¡¡ ¡¡¡ ¡¡¡ ¡¡¡ ¡¡¡ ¡¡¡ ¡¡¡ ¡¡P 1 P2 P3 P4 C1 C2 C3 C4 C5 C6
It is also relevant to use the blocks computed by our approach in order to assess the quality of the representation of the cube. More precisely, this quality is related to the following:

-the proportion of cells included in the blocks (the higher the proportion, the lower amount of cells not covered by the rules),

-the number of built blocks (the more blocks there are, the more data are heterogeneous),

-the number of blocks in comparison with the number of measure values (if several blocks are built for the same measure value, then it is not stored in a contigous manner),

-the number of overlappings between blocks and their size (the higher the number of overlapping blocks, the more mixed are the data).

Multidimensional Databases

"A data warehouse is a subject-oriented, integrated, time-variant and non-volatile collection of data in support of management's decision making process" [START_REF] Inmon | Building the Data Warehouse[END_REF]. Such systems, devoted to intensive decision-oriented querying, are different from classical relational database management systems which are not suitable in the OLAP framework. Multidimensional databases have thus been proposed in [START_REF] Codd | Providing OLAP to user-analysts: An IT mandate[END_REF]. Data are represented as multidimensional tables. No consensual definition has emerged for now concerning data representation and manipulation.

Generally, a multidimensional database is a set of hypercubes (hereafter cubes). A cube can be seen as a set of cells and a cell represents the association of a measure with one member in each dimension. Hierarchies may be defined over dimensions, for instance to describe sales in function of states and not of cities.

Definition 1 -Cube.

A k-dimensional cube, or simply a cube, is a tuple C, dom 1 , . . . , dom k , dom m , m C where -C is the name of the cube, -dom 1 , . . . , dom k are k finite sets of symbols for the members associated with dimensions 1, . . . , k respectively, -let dom mes be a finite totally ordered set of measures. Let ⊥ ∈ dom mes be a constant (to represent null values). Then

dom m = dom mes ∪ {⊥}, -m C is a mapping: dom 1 ×. . .×dom k → dom m . A cell c of a k-dimensional cube C is a (k + 1)-tuple v 1 , . . . , v k , m such that for ev- ery i = 1, . . . , k, v i is in dom i and where m = m C (v 1 , . . . , v k ). Moreover, m is called the content of c and c is called a m -cell.
Operations are defined to manipulate data cubes: selection, projection, rotation, and switch. Switch is defined to modify cube representations without altering the data. [4] details a study on cube representations. The authors propose an approach to restructure the dataset in order to obtain an appropriate representation. In this paper, we consider a fixed kdimensional cube C and a fixed representation of C, R = {rep 1 , . . . , rep k }. Given a dimension d i in C, and v 1 and v 2 in dom i , v 1 and

v 2 are said contiguous if rep i (v 1 ) and rep i (v 2 ) are consecutive integers, i.e. if |rep i (v 1 ) - rep i (v 2 )| = 1. Moreover, the interval [v 1 , v 2 ]
is the set of all contiguous values between v 1 and v 2 . In our approach, a block of C is a sub-cube of C. Note that in the case where an interval is not specified on a dimension d i , the interval is set to δ i = ALL = dom i . Therefore, we shall consider a block as defined by k intervals. Definition 4 -Block Overlapping. Two blocks overlap if they share at least one cell.

It is easy to see that two blocks

b = δ 1 × . . . × δ k and b = δ 1 × . . . × δ k overlap iff for each dimension d i , δ i ∩ δ i = ∅.
In our formalism, a slice is defined as a particular block.

Definition 5 -Slice. Let v i be a value from

dom i . A slice of C associated with v i , denoted T (v i ), is the block δ 1 × . . . × δ k such that δ i = {v i },
and for all j = i, δ j = ALL.

A slice is an hyperplan, reduced to a single row in the particular case of a twodimensional cube. Two slices defined on the same dimension d i are said to be contiguous if they are associated with two contiguous values from d i . For instance in Figure 1, the slices T (P 3) and T (P 4) are contiguous since P 3 and P 4 are contiguous in the considered representation. Cells are considered as neighbors if they share a side or a corner in the representation. For instance, in Figure 1 the cells P 2, C2, 8 and P 3, C1, 8 are neighbors.

Definition 6 -Cell Neighborhood. Two distinct cells c = v 1 , . . . , v k , m and c = v 1 , . . . , v k , m are neighbors if for every i = 1, . . . , k, |rep i (v i ) -rep i (v i )| ≤ 1.
We note that in a k-dimensional cube, a cell has 3 k -1 neighbors. Support and confidence are defined for blocks as follows. This property is used in our Apriori-like algorithm, given in the forthcoming section. Given a measure value m, a support threshold σ, a confidence threshold γ, Most-General and Most-Specific blocks with respect to m, σ and γ are defined as follows.

Definition 9 -MG-Block. Given σ a support threshold, m a measure value, and b a σ-frequent block for m, b is said to be most-general (MG) for m and the confidence threshold γ if

• conf (b, m) > γ • there does not exist a block b such that: -b is σ-frequent for m -∃j ∈ [1, k] such that δ j = ALL and δ j = ALL -∀j ∈ [1, k], j = j ⇒ δ j = δ j -conf (b , m) > γ.
Definition 10 -MS-Block. Given σ a support threshold, m a measure value, and b a σ-frequent block for m, b is said to be mostspecific (MS) for m and the confidence threshold γ if

• conf (b, m) > γ • there does not exist a block b such that: -b is σ-frequent for m -∃j ∈ [1, k] such that δ j = ALL and δ j = ALL -∀j ∈ [1, k], j = j ⇒ δ j = δ j -conf (b , m) > γ.

Block Generation

In this paper, our goal is to discover blocks b such that supp(b, m) > σ and conf (b, m) > γ, where σ and γ are respectively support and confidence thresholds specified by the user, and where m is a measure value appearing in the cube. Our method is based on a levelwise Apriori-like algorithm [START_REF] Agrawal | Mining Association Rules Between Sets of Items in Large Databases[END_REF], for scalability reasons. The proposed algorithm builds MS-blocks, which correspond to rules in which the measure value is determined by the greatest number of dimensions. However, this algorithm is also suitable for the mining of MG-blocks, which correspond to the less specific rules. The user can decide whether he (she) wants MS-or MG-blocks. The proposed method is outlined in algorithms 1 and 2. Algorithm 2 computes the blocks, provided the intervals are computed.

Regarding the computation of the intervals, we introduce the following notation: Given a measure value m and a slice T (v) associated to value v in dom i , let v be the value in

dom i such that rep -1 i (v ) = rep -1 i (v) -1. If c is a m-cell in T (v)
, then we denote by µ(c) the number of m-cells in T (v ) that are neighbors of c.

Then, the computation of the intervals for dimension d i is based on two additional thresholds denoted by µ and ν. µ is a minimum number of m-cells in T (v ) that are neighbors of a given m-cell c in T (v). ν is the minimum ratio between the number of m-cells in T (v) having at least µ m-cells neighbors in T (v ) and the total number of m-cells in T (v). The following example illustrates how µ and ν are used in our approach. In order to see how our method can actually output b 1 and b 2 , let us consider µ = 2 and ν = 60%, and let us detail the computation for the slice T (C3), assuming that T (C1) and T (C2) have already been considered. In this case, an interval of the form [C1, N IL] is under construction (N IL stands for an unknown value). At this stage of the computation, we know that C1 and C2 belong to this interval, and we have to decide for C3. To this end, we first compute supp(T (C3), 8) and, as the result is greater than 1/12, the process continues (note that, otherwise N IL would be replaced by C2). The two 8-cells in T (C3) to consider are c = P 3, C3, 8 and c = P 4, C3, 8 , and we have µ(c) = µ(c ) = 2. Since T (C3) contains two 8cells, then the ratio between the number of 8-cells in T (C3) having at least µ 8-cells neighbors in T (C2) and the number of 8-cells in T (C3) is 1, which is greater than 60%. In this case, we consider that C3 should belong to the interval under construction, and thus, its current value remains Figure 4 illustrates the construction of the blocks for value 8. Figure 5 illustrates the result of this process for all measure values [START_REF] Agrawal | Automatic Subspace Clustering of High Dimensional Data for Data Mining Applications[END_REF][START_REF] Codd | Providing OLAP to user-analysts: An IT mandate[END_REF][START_REF] Ester | Incremental Clustering for Mining in a Data Warehousing Environment[END_REF][START_REF] Inmon | Building the Data Warehouse[END_REF]75). Four blocks are built: It is also worth mentioning that the support threshold determines the minimal size of blocks while the confidence threshold determines the homogeneity of the cell values within a block. Indeed, for a given value of support σ, denoting by N the number of cells of the cube, a block can be frequent only if it contains at least σ * N cells. Moreover, for a given value of confidence γ, a block of cardinality M is kept only if it contains at least γ * M cells having the measure value m.

The algorithm is easily modified to build less specific rules. For this purpose, at each step the confidence of blocks is computed, and the scan of dimensions is stopped for each block having a sufficient confidence. In our example (Figure 5), the following sets of blocks are obtained:

B 1 = ∅, B 2 = {b 3 }, B 3 = {b 2 , b 4 }, B 4 = {b 3 }.

Rule Generation

In our model, one rule is associated with each block. When a block b j = δ 1,j × . . . × δ k,j does not overlap any other block (B j = ∅), the generated rule for b j is as follows:

If d 1 ∈ δ 1,j and . . . and d k ∈ δ k,j Then m j where m j is the associated measure value. For instance, the rule built for the block b 1 (Figure 1) is:

Algorithm 1: Computation of L 1 (m) L 1 (m) ← ∅ ; foreach dimension d i i=1,. . . ,k do intervals(m, i) ← ∅ ; currentInterval ← [N IL, N IL] ; foreach j = 1, . . . , |dom i | do compute supp(T (rep -1 i (j)), m) ; if supp(T (rep -1 i (j)), m) ≤ σ then if currentInterval = [α
, i) ← intervals(m, i) ∪ {[rep -1 i (α), rep -1 i (j -1)]} ; currentInterval ← [j, N IL] ; L 1 (m) ← L 1 (m) ∪ intervals(m, i) ; B 1 (m) ← {δ 1 × . . . × δ k | (∃i, δ i ∈ intervals(i, m)) and (∀j = i, δ j = ALL)}
If CITY ∈ [C1, C1] and PRODUCT ∈ [P 1, P 2] Then 6
If an overlapping is detected, then a fuzzy rule is built. Membership functions can be defined in different manners. In our approach, we define the membership function based on the number of overlapping blocks: for dimension values where there is no overlapping, the degree is 1, and if there are p overlapping blocks, then the degree is1 p . Algorithm 3 describes this computation.

Γ P RODU CT,3 (P ) =    1 if P ∈ [P 1, P 2] 1 3 if P ∈ [P 3, P 3] 0 otherwise Γ CIT Y,3 (C) =    1 2 if C ∈ [C3, C3] 1 2 if C ∈ [C4, C4] 0 otherwise
Note that the intervals [C3, C3] and [C4, C4] can be merged into the single one [C3, C4]. Figure 6 illustrates the membership functions built for each block (or rule).

i,j = {[α 1 i,j , β 1 i,j ], . . . , [α r i,j , β r i,j
]} be the set of all sub-intervals of [α i,j , β i,j ] ; Associate membership function Γ i,j to [α i,j , β i,j ] as defined by: foreach sub-interval [α s i,j , β s i,j ] ∈ Int i,j do foreach value v from dom i do Γ i,j (v) = 0 if v ∈ [α s i,j , β s i,j ]

the goals are the same, it is not possible to apply such methods due to problems of scalability and because also of the multidimensional nature of data.

Segmentation methods (e.g. clustering) have been proposed in the multidimensional context [START_REF] Agrawal | Automatic Subspace Clustering of High Dimensional Data for Data Mining Applications[END_REF][START_REF] Ester | Incremental Clustering for Mining in a Data Warehousing Environment[END_REF]. In [START_REF] Gyenesei | Multidimensional Partitioning of Attribute Ranges for Mining Frequent Fuzzy Patterns[END_REF], the authors study the generation of fuzzy partitions over numerical dimensions. However, these propositions are not related to the measure value and thus are different from our work where the measure value is the central criterion.

[9] aims at compressing data cubes. However there is no consideration on cube representations and homogeneous blocks generation.

[3] proposes a method to divide cubes into regions and to represent those regions. However, the authors aim at representing the whole cube. They use statistical methods to construct an approximation of the cube, while we aim at discovering relevant rules as fast as possible which may not cover the whole cube.

In [START_REF] Wang | Condensed Cube: An Effective Approach to Reducing Data Cube Size[END_REF], the authors propose the concept of condensed data cube. However, the authors aim at condensing the cube without loss of information while we aim at displaying relevant information to the user, which may be a partial representation of data.

Conclusion

In this paper, an efficient method for summarizing and visualizing multidimensional data is proposed. In our approach, blocks of homogeneous data are built and the block overlapping is managed using fuzzy sets. Each block is associated with one rule (fuzzy or not). The set of these rules is meant to summarize the cube. Further work include the implementation of our method, the introduction of fuzziness within blocks in order to discover blocks of almost the same value instead of homogeneous blocks, and the study of the pre-treatments in order to organize cubes in a relevant manner before mining blocks.

Figure 1 :

 1 Figure 1: Data cube and associated blocks Building such blocks facilitates visualization.The more the relevance of the blocks, the more the representation quality. In[START_REF] Choong | Computing Appropriate Representation for Multidimensional Data[END_REF], the authors study different manners to represent data cubes. In particular, the authors highlight the fact that some representations are more relevant than others since they allow to display the cube according to user-specified criteria. In this work, we consider as relevant the configurations where same measure values are grouped. Another organization is possible, as described in[START_REF] Choong | Computing Appropriate Representation for Multidimensional Data[END_REF], where the data are organized so that the measure is ordered in an increasing manner over all the dimensions. However, it is an NP-hard problem to automatically build such a representation of a cube. In this paper, the organization of data is not altered before the block identification, although it could be interesting.

Figure 2 :

 2 Figure 2: Inverting dimension values.

Definition 2 -

 2 Representation. A representation of a cube C is a set R = {rep 1 , . . . , rep k } where for every i = 1, . . . , k, rep i is a one-to-one mapping from dom i to {1, . . . |dom i |}.

Figures 1 and 2

 2 Figures 1 and 2 display two different representations of the same cube.

Definition 3 -

 3 Block. A block b is a set of cells defined over a cube C by b = δ 1 × . . . × δ k where δ i are intervals of contiguous values from dom i , for i = 1, . . . , k.

Definition 7 -

 7 Support. The support of a block b from C for a measure value m is defined as: supp(b, m) = # cells containing m in b # cells in C . Considering a user-given minimum support threshold σ and a measure value m, a block b such that supp(b, m) > σ is called σfrequent for m. Note that the support is antimonotonic, meaning that for all blocks b, b and for each measure value m: b ⊆ b ⇒ supp(b, m) ≤ supp(b , m).

Definition 8 -

 8 Confidence. The confidence of a block b for a measure value m is defined as: conf (b, m) = # cells containing m in b # cells in b .

Figure 3 :

 3 Figure 3: The Cube for Example 1

Example 1

 1 Consider the cube C shown in Figure 3 in which blocks are to be computed according the support and confidence thresholds equal to 1/12 and 50% respectively. Let b be the block [P 1, P 4]×[C1, C5]. Clearly, supp(b, 8) > 1/12, but conf (b, 8) ≤ 50%, showing that b cannot be output by our method. On the other hand, for the measure value 8, the blocks b 1 = [P 1, P 2] × [C4, C5] and b 2 = [P 3, P 4] × [C1, C3] have their supports and confidences greater than 1/12 and than 50% respectively.

Figure 4 :

 4 Figure 4: Occurences of measure value 8

b 1 =

 1 [P 1, P 2] × [C1, C1] for value 6 b 2 = [P 3, P 4] × [C1, C3] for value 8 b 3 = [P 1, P 3] × [C3, C4] for value 5 b 4 = [P 3, P 4] × [C4, C6] for value 2 Note that there are two overlappings (between b 2 and b 3 , and between b 3 and b 4 ).

Figure 5 :

 5 Figure 5: Output of Algorithm 1 As seen previously, blocks output by the above algorithm may overlap. In our approach, overlappings are computed based on the fact that two blocks overlap if and only if all intervals defining these blocks have a nonempty intersection. The proposed algorithm, which is omitted here, associates each block b j in B with the set B j of all blocks b in B such that b = b j and b overlaps b j .

  , N IL] then /*close the current interval at position j -1 and set the current interval to the empty interval */ intervals(m, i)← intervals(m, i) ∪ {[rep -1 i (α), rep -1 i (j -1)]} ; currentInterval ← [N IL, N IL] ; else if currentInterval = [NIL, N IL] then /*start a new current interval at position j */ currentInterval ← [j, N IL] ; else /* currentInterval = [α, N IL] */ ; Let neighbors(j, m) be the number of cells in T (rep -1 i (j)) having at least µ neighbors in T (rep -1 i (j -1)) containing m ; Let M be the number of cells in T (rep -1 i (j)) containing m ; if (neighbors(j, m)/M ) ≤ ν then /*close the current interval at position j -1 and start a new current interval at position j */ intervals(m

Algorithm 2 :

 2 Discovery of MS-blocksData : data cube C, σ: support threshold, γ: confidence threshold, ν: homogeneity threshold. Result : set of blocks B associated with C foreach measure value m from C do Compute L 1 (m) ;for l = 2 to k do B l (m) ← ∅ ; Generate from L i l-1 candidates δ i1 × . . . × δ i l such that ∀p, p ∈ [1, l], i p = i p ; Let L l (m) be this set ; Pruning: Delete from L l (m) all candidates δ i1 × . . . × δ i l such that there exists p ∈ {1, . . . , l} such that δ i1 × . . . × δ ip-1 × δ ip+1 × . . . × δ i l is not frequent. ; foreach remaining candidate δ i1 × . . . × δ i l doLet b be the block δ 1 × . . . × δ k where δ p = δ pj if dimension d p has been treated andδ p = ALL otherwise. if supp(b, m) ≤ σ then remove δ i1 × . . . × δ i l from L l (m) else B l (m) ← B l (m) ∪ {b}Let B(m) be the set of all MS blocks b in l=k l=1 B l (m) such that conf (b, m) > γ B ← m B(m) if the data cube must be visualized then foreach b ∈ B do replace the values of all cells of b by the measure value associated to b. /*Cells not belonging to a block are set to null. (Note that a cell may contain several values.)*/

  For instance, consider block b 3 = [P 1, P 3] × [C3, C4] with B 3 = {b 2 , b 4 }. For dimension P RODU CT, the interval δ P RODU CT,3 = [P 1, P 3] is divided into two sub-intervals [P 1, P 2] and [P 3, P 3]. Applying the same method on dimension CIT Y , we obtain:

Figure 6 :

 6 Figure 6: Fuzzy Set Construction

  Computation of Fuzzy SetsData : n sets of blocks B j Result : Fuzzy Sets foreach set of blocks B j = {b j1 , . . . , b j l } with b jp = δ 1,jp × . . . × δ k,jp and δ r,jp = [α r,jp , β r,jp ] (r = 1, . . . , k) do foreach dimension d i do foreach q = 1, . . . , l do α i,jq ← max(α i,j , α i,jq ) ; β i,jq ← min(β i,j , β i,jq ) ; Sort the set {α i,jq | q = 1, . . . , l} ∪ {β i,jq | q = 1, . . . , l} ; let Int

	Algorithm 3:
	5 Related Work
	Research work on (fuzzy) image segmentation
	may appear as related works [10]. Although

(1+ns) otherwise n s : # blocks b jp ∈ B j s.t. [α s i,j , β s i,j ] ⊆ [α i,jp , β i,jp ]

Acknowledgements

The authors acknowledge the French Embassy in Malaysia for its support.