A Carcaterra

New concepts in damping generation and control: theoretical formulation and industrial applications

Introduction

These notes are finalized to a particular study of the damping mechanism in Hamiltonian systems, characterized indeed by the absence of any energy dissipation effect. It is important to make a clear distinction between the two previous concepts, since they seem to be somehow contradictory. A Hamiltonian system is characterized by an invariant total energy (the Hamiltonian H) that is equivalent to state any energy dissipation process is absent. This circumstance, especially from an engineering point of view, leads to the wrong expectation that the motion of any part of such a dissipation-free system, subjected to some initial conditions, maintains a sort of constant amplitude response. This is, although unexpectedly, a wrong prediction and the "mechanical intuition" leads, in this case, to a false belief. It is indeed true the converse: even in the absence of any energy dissipation mechanisms, mechanical systems can exhibit damping, i.e. a decay amplitude motion.

This statement makes clear how the two concepts rely on completely different properties of the Hamiltonian systems.

This fact is the root of very fundamental physical properties of mechanical systems, and touches charming and thorny questions about them. It is a matter of fact that looking at the atomic scale of mechanical systems, they can be described through the use of Hamiltonian equations: an atomic lattice vibrates and moves energy along its structure without any dissipation mechanisms. But it is a trivial consideration about every day physical world, that any energy released at macro-scale to an atomic lattice tends to disappear at that scale, to be moved into the micro-scale vibrations of the atoms. Damping is observed. This process is the base for converting mechanical energy (large scale motion) into heat (small scale thermal vibration), and without any energy loss. The natural tendency to produce this process is at the base of the macroscopic irreversibility and at the root of the second principle of thermodynamics itself.

As it appears, the questions behind the distinction between damping and dissipation involves fundamental aspects of physical systems, that have been the subject of diatribes, debates and fascinating investigations starting from the end of the nineteenth century up to the present days, in the field of theoretical physics as well as in mathematical physics. Famous names are involved as those of Boltzmann, Loschmidt, Curvelwell, the Ehrenfests, Poincarè, Zermelo, Prigogine and many others, turning around the Boltzmann's H-theorem.

The investigation contained in the present notes does not attack directly the previous questions (wisely), but a problem somehow close to the previous, and the related possibilities for engineering applications(!). As a consequence, the obtained results seem to be theoretically intriguing (aesthetic in science is a luxury relatively down market) and practically usable (engineers are sensitive). Optimists would say a good compromise.

But, aside these considerations, let us illustrate the main point we have here.

More precisely, we consider a partition of the Hamiltonian system S into two subsystems: we select, among its N degrees of freedom, one of them x, indicated as master. The remaining part of the system, consisting of N -1 degrees of freedom x, is called the hidden or unmonitored part of the system. With this view of the problem, we are interested in the following analysis.

Let the Lagrangian function L of S be:

L (x, ẋ) = L x (x, ẋ) + L mix x, ẋ, x, ẋ + L hid x, ẋ
The motion of the whole system S is governed by the minimum principle:

δ t2 t1 L (x, ẋ) dt = 0
to which correspond the Euler-Lagrange equations:

d d t ∂L (x, ẋ ) ∂ ẋ - ∂L (x, ẋ) ∂x = 0 (1)
providing the solutions x * , x * .

The motion x * would be here derived directly from the solution of a reduced problem, by means of a new equation of motion written in terms of x only in the form:

d d t ∂ L x (x, ẋ) ∂ ẋ - ∂ L x (x, ẋ) ∂ x + Q (x, ẋ, ẍ) = 0 (2)
where the form for Q depends on the structure of the Lagrangian L.

Comparing the complete set of equations for S and that for x, the following replacement is operated:

V x, ẋ, x, ẋ = d d t ∂L mix x, ẋ, x, ẋ ∂ ẋ - ∂L mix x, ẋ, x, ẋ ∂x replaced by =⇒ Q (x, ẋ, ẍ)
Thus, we desire to determine, if possible, the form of Q such that equations ( 1) and ( 2) provide the same solution x* for x. It is shown in the sections ahead that in general it is not strictly possible but, under some restrictive hypotheses, the difference between Q and V can be small, and negligible, at least in a prescribed time window. The nature of the term Q under such hypotheses, can be also nicely interpreted as a dissipative effect for the system x; it means:

L x (x, ẋ) = T x (x, ẋ) -U x (x) ⇒ Ḣx = d dt [T x (x, ẋ) + U x (x)] < 0 (3)
Therefore, such analysis explains well the chance of observing damping in the motion of x even if dissipation in the system S is absent i.e.: Ḣ = d dt [T (x, ẋ) + U (x)] = 0 [START_REF] Weaver | The effect of an undamped finite degree of freedom 'fuzzy' substructure: numerical solution and theoretical discussion[END_REF] being clear how equations [START_REF] Maidanik | Induced damping by a nearly continuous distribution of a nearly undamped oscillators: linear analysis[END_REF] and ( 4)are not contradictory. The question of the form of Q is approached in Section 2, namely subsections 1, 2 and 3, where it is illustrated how Q can include a dissipation term. Subsection 4 studies the time window within which the substitution of V by Q is permitted because of the small error.

Subsection 5 approach an inverse problem through a variational technique: design the system S such that Q generates the fastest energy transfer from x to the unmonitored system. Subsection 6 studies examples of application of the theory presented in 5.

Finally, last but not least, Section 3 describes an industrial application of the theory outlined, showing how the theoretical speculation can dress the clothes of a design procedure for the production of a real innovative category of damping devices.

Therefore, the main point we have here is that very effective damping properties can be indeed obtained by coupling a main structure to a secondary structure that serves as energy storage. This opportunity opens a new way to damping control: damping is produced not acting on a local dissipation process, a weak controllable phenomenon, but it can be shaped by the designer modifying the purely elastic response of a secondary part of the structure. This result can be achieved by purely mechanical means or by electromechanical devices. The theory developed here applies to both, and in both cases the ability in controlling the damping properties of the primary structure does not rely on any local dissipation mechanism.

More precisely, it has been shown how the motion of a principal structure, called master, can be damped as an effect amounting to an energy transfer process by which the master energy is moved to a set of resonators attached to it (1; 2; 3; 4; 5; 6; 7; 8). The phenomenon, sometime called apparent damping, is intriguing considering that this spontaneous energy transfer can have irreversible characteristics, i.e. the energy can be permanently transferred from the master to the set of resonators. This result is always observed when the number of degrees of freedom within the set is infinite, or actually very large (9; 10; 11; 20; 21; 22; 23; 26; 27; 28) and not for small number of oscillators, where recursive phenomena are observed (8; 12). In some cases, the same irreversible energy transfer can be indeed predicted considering statistical ensemble average over a population of similar structures [START_REF] Carcaterra | Ensemble energy average and energy flow relationships for nonstationary vibrating systems[END_REF]. Moreover, a quasi-irreversible energy transfer can be produced even with a finite number of resonators within the attachment, when selecting special distribution of the natural frequencies of the attached resonators, or alternatively, introducing non-linear or parametric effects within the set (14; 15; 16; 17; 27; 28). In some cases the attachment can have an electrical nature (24; 25). In (38; 39; 40; 41; 42; 43; 44; 45; 46; 47) it is shown how an electrical passive network, a sort of electrical-double of the master mechanical structure, can absorb very effectively the energy pumped in it. An energy transfer from the mechanical to the electrical part of the system initially takes place; later, energy is actually dissipated by resistive effects into the electrical energy storage. Also in this case, the observed damping is not a direct consequence of a local dissipation, but rather it is due to a fast energy displacement form the master part of the system to a sacrificial attachment.

Further interesting readings on the subject can be found in (18; 19; 29; 30; 31; 32; 33; 34; 35; 36; 37).

Pseudo-Dissipative Systems: an Outline of the Theory

This section outlines an original theory of pseudo-dissipative structures, Ndimensional Hamiltonian systems for which a certain degree of freedom x, called master, is separated from the remaining N -1 dofs, the hidden variables, and studied apart from them. The effect of the interaction between x and the hidden dofs is taken into account by terms dependent only on x, following the procedure outlined in these notes. Namely the attention is focused on three main aspects: (i) collapse this interaction effect into simple terms, (ii) provide for them a physical interpretation, (iii) discuss in depth the limit of these simplifications.

The analysis starts with static systems. For them the effect of the hidden variables is not a pseudo-dissipation but rather an additional elastic restoring force depending on the displacement of x and on the system of forces acting on the hidden variables. However, the static system allows to understand simply the way the hidden variables play their effect in the equation of motion of the master variable, especially looking at the Hamiltonian formulation of the reduced problem. This presents a formal analogy with the dynamic case, the one of main interest for the theory presented in this course.

The Hamiltonian presentation of the results permits to understand how the terms depending on the hidden variables and on the mixed terms, are indeed replaced by simpler equivalent terms in the new Lagrangian function of the reduced system. On the other hand, the variational approach plays a key role in determining the conditions for a permanent energy storage within the master variable x.

Static systems: the effect of the hidden variables

A static system S is described through the set of variables

x ≡ x 0 , x 1 , . . . , x N x = {x 0 , x 1 , . . . , x N } T
Assume S is controlled by the set of equations:

∂L (x) ∂x k = 0, k= 0, 1, 2, . . . , N
the Euler-Lagrange equations related to the variational principle:

δL (x) = 0
where the characteristic lagrangian function L for S is dependent on x in a quadratic fashion

L (x) = 1 2 x T K x -f T x, K =   k 00 . . . k 0N . . . . . . . . . k N 0 . . . k NN   = k 00 kT 0 k0 K K =   k 11 . . . k 1N . . . . . . . . . k N 1 . . . k NN   =   kT 1 . . . kT N   , k0 =    k 01 . . . k 0N    =    k 10 . . . k N 0   
Accordingly to our assumption the equations for S are:

∂ L (x) ∂ x k = 0, → K x = f
These can be written separating the degree of freedom x, the master, from the others as:

         k 00 x + kT 0 x = f 0 k 10 x + kT 1 x = f 1 . . . k N 0 x + kT N x = f N or    k 00 x + kT 0 x = f 0 k0 x + K x = f
where:

x =    x 1 . . . x N    , f =    f 1 . . . f N   
We can eliminate from the equation for x the set of other variables (hidden) derived as the solution of the matrix equation, obtaining a single equation for x, where the hidden variables disappear:

   k 00 x -kT 0 K -1 k0 x = f 0 -kT 0 K -1 f x = K -1 f -k0 x
The determined equation for x deserves some comments. Consider explicitly the lagrangian of the original problem, decoupling the dependency of it upon x and the other variables:

L (x) = 1 2 x T Kx -f T x = 1 2 k 00 x 2 + 1 2 x kT 0 x + 1 2 x kT 0 x + 1 2 xT Kx -f 0 x -f T x
We can distinguish three kind of terms: terms directly dependent on x, terms dependent on the hidden variables, and mixed terms:

L (x) = L X (x) + L mix (x, x) + L hid (x) L X (x) = 1 2 k 00 x 2 -f 0 x L mix (x, x) = x kT 0 x L hid (x) = 1 2 xT K x -f T x
It comes out the equation for x can be determined by modifying the lagrangian function L of the complete system retaining only the direct terms and substituting the mixed terms and the terms depending on the hidden variables by an interaction potential D(x) and an external force potential N x:

L X (x) = L X (x) + D (x) + N x where: D (x) = 1 2 kT 0 K -1 k0 x 2 N = kT 0 K -1 f

Dynamic systems

The Hamiltonian system S is described through the set of variables x (t) ≡ x 0 (t) , x 1 (t) , . . . , x N (t)

x (t) = {x 0 (t) , x 1 (t) , . . . , x N (t)} T and governed by the set of equations:

d d t ∂ L (x, ẋ ) ∂ ẋk - ∂ L (x, ẋ) ∂ x k = 0, k= 0, 1, 2, . . . , N
the Euler-Lagrange equations related to the well-known variational Hamilton principle:

δ t2 t1 L (x, ẋ) dt = 0
where L is the quadratic form:

L (x, ẋ) = 1 2 ẋT M ẋ + 1 2 x T Kx -f T x K =   k 00 . . . k 0N . . . . . . . . . k N 0 . . . k NN   = k 00 kT 0 k0 K K =    k 11 . . . k 1N . . . . . . . . . k N 1 . . . k NN    =    kT 1 . . . kT N    , k0 =    k 01 . . . k 0N    =    k 10 . . . k N 0    , M =    m 0 0 . . . 0 0 m 1 . . . 0 0 0 . . . m N   
This partitioning of matrices is useful ahead. Accordingly to our assumption, the equations for S are:

d d t ∂ L (x, ẋ ) ∂ ẋk - ∂ L (x, ẋ) ∂ x k = 0, → Mẍ + Kx = f
We can separate x from the other degrees of freedom as:

           m 0 ẍ + k 00 x + kT 0 x = f 0 m 1 ẍ1 + k 10 x + kT 1 x = f 1 . . . m N ẍN + k N 0 x + kT N x = f N or    m 0 ẍ + k 00 x + kT 0 x = f 0 M ẍ + k0 x + K x = f
where:

x =    x 1 . . . x N    , M =   m 1 0 . . . 0 0 m 2 . . . 0 0 0 . . . m N   , f =    f 1 . . . f N   
The lagrangian function can be written separating the contribution of x from the others as:

L (x, ẋ) = 1 2 ẋT M ẋ + 1 2 x T Kx -f T x = = 1 2 m 0 ẋ2 + 1 2 k 00 x 2 + 1 2 x kT 0 x + 1 2 x kT 0 x + 1 2 ẋT M ẋ+ + 1 2 xT Kx -f 0 x -f T x
The Laplace domain form of the equations of motion reads (capital symbols for transformed quantities):

           m 0 s 2 X + k 00 X + kT 0 X = F 0 m 1 s 2 X 1 + k 10 X + kT 1 X = F 1 . . . m N s 2 X N + k N 0 X + kT N X = F N or    m 0 s 2 X + k 00 X + kT 0 X = F 0 s 2 M X + k0 X + K X = F
As for the static case, eliminating the hidden variables from the equation for x yields:

     m 0 s 2 X + k 00 X -kT 0 s 2 M + K -1 k0 X = F 0 -kT 0 s 2 M + K -1 F X = s 2 M + K -1 f -k0 X
Transforming the first equation back to the time domain, one obtains the three equivalent forms:

m 0 ẍ(t) + k 00 x(t)- 1 2π ∞ 0 kT 0 G(s) k0 X(s)e st dt = f 0 - 1 2π ∞ 0 kT 0 G(s) F(s)e st dt m 0 ẍ(t) + k 00 x(t) + ∞ -∞ kT 0 g(t -τ ) k0 x(τ )dτ = f 0 + ∞ -∞ kT 0 g(t -τ ) f (τ )dτ m 0 ẍ(t) + k 00 x(t) + kT 0 g(t) k0 * x(t) = f 0 + kT 0 g(t) * f (t)
where

g(t) = - 1 2π ∞ 0 s 2 M + K -1 e st dt ; G(s) = s 2 M + K -1
These equations produce a clear qualitative picture of the nature of the interaction effect between x and the hidden variables. More precisely, note that: 1. The reduced equation of the motion for x becomes an integral-differential equation, meaning the interaction with the remaining part of the system, the hidden part, amount to a memory effect (integral term); this is the dynamic counterpart of the terms related to the interaction potential D(x) appearing in the static case;

2. Additionally, the forces applied to the hidden part of the system appears as a known forcing term, that amounts to a noise effect; this is the dynamic counterpart of the static external force potential Nx ; 3. The kernel G apparently is a rational function in terms of s; namely it can be expressed as a ratio of two polynomials P and Q (the degree of P being larger than that of Q, because of the causality principle):

[G(s)] ij = s 2 M + K -1 ij = P ij (s) Q(s)
4. The kernel G, since the system is Hamiltonian, is a real function, and does not contain any imaginary part, as an effect of the absence of any real dissipation in the motion of the system. 5. In the previous analysis, using the Laplace transform, it has been tacitly assumed that all initial conditions are set to zero (generalization is relatively easy). As the simplest example, a nice exercise it to apply the previous analysis to a two degrees of freedom system; one can easily verify that the previously discussed equations becomes:

m 0 ẍ(t) + (k 1 + k 2 )x(t) + k 2 2 g(t) * x(t) = f 0 = k 2 g(t) * f 1 (t) M = m 0 0 0 m 1 , K = k 1 + k 2 -k 2 -k 2 k 2 G(s) = 1 m 1 s 2 + k 1 , g(t) = H(t) sin k 2 m 2 t 2.
3 Pseudo-dissipative effects: motion about ω 0

The previous analysis shows how following the motion of only one selected degree of freedom of S, its motion becomes controlled by an integraldifferential equation, the integral part amounting to the interaction of x with the unmonitored part of the system. However, a more clear and physically interesting picture of the problem arises from a more detailed analysis of this integral terms. If some approximations are made, it discloses sharp characteristics about the nature of the forces they represent.

To enlightening these properties, our following analysis offers two distinct mathematical approaches, that leads to look at different perspectives for the interaction terms. Both the approaches are based on an approximation about the integral term involving G. The idea is to consider systems having a dominant frequency in their response, say ω 0 . This happens when the natural frequencies of S are all located in a rather narrow frequency bandwidth. In physical systems this case is actually met in many cases of interest and it is also a good approximation to approach the energy sharing process in the general case. For engineering system this can indeed intended as a design configuration purposely obtained to determine a desired energy sharing effect.

However, both of the proposed techniques permit to extract information about the energy exchange between x and the hidden variables. More precisely the first approach passes through a Taylor expansion of the kernel G in terms of s, while the second uses a less intuitive strategy, but more powerful, we named integral-Padè expansion of the kernel.

In the next sections the two mentioned approaches are developed in detail.

Pseudo-dissipative effects: Taylor expansion of the kernel G To carry on our point of view about unmonitored dynamic systems, it helps introduce some additional considerations and then some additional hypotheses on our system.

The kernel G can be written using a Taylor expansion in the complex plane in terms of s, valid within a circle around the complex s 0 :

G (s) = G (s 0 ) + G (s 0 ) (s -s 0 ) + R (s) R (s) = (s -s 0 ) N N ! d (N ) G ds N s=ζ = 1 2πi γ(s0) G(s) (s -s 0 ) N +1 ds
Note that this expansion of G produces, neglecting the remainder, a linear term G (s 0 )s plus a constant term G(s 0 ) -G (s 0 )s 0 . In general G(s 0 ), G (s 0 ) are complex quantities and it is useful to analyse them in depth. We assume:

s 0 = σ 0 + jω 0 , s= jω
This means: the center of the series expansion is in general on the complex plane, while the series is evaluated only along the imaginary axis. Since the expansion is made in the frequency domain (ω) about ω 0 , we call the investigated response "motion about ω 0 ". The physical sense of this assumption and several cases of applications will be the subject of following subsections.

Neglecting the remainder one gets the equation:

G ≈ Re(G 0 -σ 0 G 0 )-(ω-ω 0 )Im(G 0 )+j [Im(G 0 -σ 0 G 0 )+(ω -ω 0 )Re(G 0 )]
This expression leads, in general, to four different force contributions in the equation of motion, that can be categorized as: Dissipative actions: (i) viscous damping: jωRe (G 0 ) (ii) hysteretic damping:

j [Im (G 0 -σ 0 G 0 ) -ω 0 Re (G 0 )] Conservative actions: (iii) elastic: Re (G 0 -σ 0 G 0 ) + ω 0 Im (G 0 ) (iv) gyroscopic effect: -ωIm (G 0 )
As it is clear two of them are dissipative, the other two conservative. The last term is here called gyroscopic because this is a velocity dependent conservative contribution, in that the net energy added to the system along one period is zero, and it adds and subtracts the same amount of energy in half of the period, that is the characteristic of gyroscopic effects. For example, in the Euler equations for the rigid body motion, the components of the angular velocity are conservative terms coupling the differential equations in terms of angular velocity components derivatives.

The hysteretic damping is indeed a velocity independent effect frequently met (at least in the frequency domain) in structural dynamics, to represent inherent dissipation of the material.

To put the equation for x into a more concise form, let:

c eq = -kT 0 Re(G 0 ) k0 η eq = - kT 0 [Im(G 0 -σ 0 G 0 ) -ω 0 Re(G 0 )] k 00 k eq = -kT 0 Re(G 0 -σ 0 G 0 + ω 0 Im(G 0 )) k0 g eq = -kT 0 Im(G 0 ) k0 F + (•) = ∞ 0 (•)e jωt dt one sided Fourier tr.
Considering again the Laplace domain equation for x:

m 0 s 2 X + k 00 X -kT 0 G k0 X = F 0 -kT 0 G F Assume for
the sake of simplicity the forces on the hidden variables zero (zero noise); with the previous determined expression for G, and with s = jω it becomes:

-m 0 ω 2 + (k 00 + k eq ) + jωc eq + jk 00 η eq + ωg eq X(jω) = F 0 It is useful to treat the gyroscopic term in a slightly different manner, as:

-m 0 ω 2 + (k 00 + k eq ) + jωc eq + jk 00 η eq + j(jωg eq ) X(jω) = F 0 One can notice that, under the assumptions made, the interaction between x and the others degrees of freedom generates: (i) pseudo-damping viscous and hysteretic effects, that amount to an energy release (or energy absorption, depending on the algebraic sign) from the coordinate x to the hidden variables, (ii) an additional restoring elastic force, (iii) a gyroscopic term.

Let now go back to the time domain representation of this equation. Introduce the relationships:

x (t) = 1 2 π +∞ -∞ X (jω) e jωt dω = F -1 (X) x (t) = 1 2π +∞ -∞ -j sign (ω) X (j ω) e jωt dω x + (t) = 1 2π +∞ 0 X (jω) e j ωt dω = F -1 + (X) = 1 2 [x (t) + j x (t)]
where the last two equations stands for the so called analytic signal and Hilbert transform of x, respectively. If the one sided inverse Fourier transform is applied to the equation of motion, one has: m 0 ẍ+ + (k 00 + k eq ) x + + c eq ẋ+ + jk 00 η eq x +jg eq ẋ+ = f 0+

Taking for this complex valued equation only the real part, one finally obtains: m 0 ẍ + (k 00 + k eq ) x + c eq ẋ + k 00 η eq x + g eq ẋ = f 0 (5) Equation ( 5) is one of the central result of the present theory. This shows that, under some hypotheses related to the remainder term, a simple and strict form for the interaction between the variable x and the hidden part of the system is born. A nice distinction can be operated between reactive forces, elastic and gyroscopic, and active forces, viscous and hysteretic. Note that among the terms depending on x, the one involving the derivative is related to dissipation, while the opposite is true for those depending on the Hilbert transform.

This result is derived under the hypothesis of keeping only the first order term in the Taylor expansion for G. This can be a reasonable approximation when the system response presents a rather narrowband spectrum concentrated around a certain reference frequency ω 0 .

Under a substantially equivalent hypothesis, a different approach can be developed that, although less intuitive, leads to a better development of our theory.

Pseudo-dissipative effects: integral-Padè expansion of the kernel G It is interesting to start a different approach to the analysis of G. The one presented in the previous section, based on the Taylor expansion, can nicely represent the interaction terms in a suggestive physical manner, but the analysis of the remainder R(s) in that form, does not allow a straightforward prediction of the effects produced when neglecting it, especially in the time domain. In other words we can expect that x is not controlled strictly by equation ( 1), and it can depart from it as an effect of the presence of R(s).

The present approach, based on a completely different idea, produces similar but simpler results in terms of the physical picture of the interaction term, but, additionally, it provides a new point of view about the pseudodissipative effect based on the analysis of the remainder term.

More precisely, we can put a list of key questions about the remainder: Q1. What is the effect of the remainder on the motion of x? Q2. Can the time domain counterpart of the remainder disclose some elements answering Q1? Q3. In which manner the form of the remainder depends on the system physical properties? Q4. The presence of the pseudo-dissipative terms, i.e. the energy released from x to the hidden variables, is a permanent property of the interaction or the existence of a not negligible remainder confines this effect in some restricted time window? Q5. There exist physical systems, purposely selected or designed, for which the effect of the remainder can be controlled, for example making it as small as possible? Q6. There is any chance to use answer to Q5 to design special engineering devices? The compiled list opens probably the core question of our problem. In fact the results we found about pseudo-dissipation have their fundament only when the effect of the remainder term is, in some sense known as unimportant. Moreover, there are two additional points that can produce a breakthrough in our theory: (i) the possible existence of particular systems exhibiting permanent pseudo-dissipative effects, as those described for equation [START_REF] Pierce | Fundamental structuralacoustic idealization for structure with fuzzy internals[END_REF], pushes our analysis toward a new border: that of the irreversible process in Hamiltonian systems, a subject of great charm as it has been outlined in the initial survey; at the same time, and not surprisingly, this circumstance drives us towards engineering design of systems able to pilot their energy to enhance their dynamic performances.

For these reasons we go more in depth in analyzing the nature of the remainder term. Starting from the Fourier domain equation for x:

-m 0 ω 2 X + k 00 X -kT 0 G (ω) k0 X = F 0 and using the Padè expansion for G (a real valued function of ω):

G ij = 2 N n=1 ϑ (i,j) n jω -jω n , G = 2N n=1 Θ n jω -jω n kT 0 G k0 = 2N n=1 kT 0 Θ n k0 jω -jω n = N n=1 α n jω -jω n + α * n jω + jω n α n = lim ω→ωn j (ω -ω n ) kT 0 G (ω) k0 = jβ n
where the poles of the transfer function occur in imaginary conjugate pairs because of the absence of any real dissipative effect in the (Hamiltonian) system.

The equation of motion becomes:

-m 0 ω 2 + k 00 X (ω) -X (ω) N n=1 α n jω -jω n + α * n jω + jω n = F 0
Taking the inverse Fourier transform and introducing the Heavisdie function H(t):

m 0 ẍ (t) + k 00 x (t) -x (t) * H (t) N n=1 α n e jωnt + α * n e -jω n t = F 0 m 0 ẍ (t) + k 00 x (t) -2x (t) * H (t) N n=1 Re α n e jωnt + α * n e -jω n t = F 0 m 0 ẍ (t) + k 00 x (t) + 2x (t) * H (t) N n=1
β n sin ω n t = F 0 Let:

S(t) = 1 N N n=1 β n sin ω n t
With this positions, the equation for x becomes:

m 0 ẍ (t) + k 00 x (t) + 2N x * [HS] = F 0 m 0 ẍ (t) + k 00 x (t) + 2N ∞ 0 x (τ ) H (t -τ ) S (t -τ ) dτ = F 0
The nature of the interaction between x and the hidden variables expressed through x * [HS], is difficult to handle because of the discrete summation appearing in S. An interesting chance comes indeed from the substitution of S by a suitable integral, that makes the physical interpretation of the interaction term easier. In fact, the core of our analysis assumes the system's eigenfrequencies are "dense enough" over the frequency axis; this means that all the natural frequencies belong to the same narrow frequency bandwidth. Under this condition it is reasonable to replace the summation over n of the Padè expansion by an integral I(t):

I (t) = 1 0 β (ξ) sin Ω (ξ) t dξ S(t) = N n=1 β n sin ω n t∆ξ ∆ξ = 1 N I(t) ≈ S(t)
where ξ is a dummy variable that belongs to the interval [0,1]. Introducing the remainder r(t):

I(t) = S(t) + r (t)
the equation of motion takes the form:

m 0 ẍ (t) + k 00 x (t) + 2N x * [HI] + x * [Hr] = F 0 m 0 ẍ (t) + k 00 x (t) + 2N ∞ 0 x (τ ) H (t -τ ) I (t -τ ) dτ + + 2N ∞ 0 x (τ ) H (t -τ ) r (t -τ ) dτ = F 0
meaning the equation is rewritten as:

m 0 ẍ (t) + k 00 x (t) + 2N x * [HI] + ε (t) = F 0 (6)
Equation ( 2) is the new starting station for our analysis. A physical interpretation of the term [HI] should be provided, together with an estimate for ε that is the force generated by the remainder term.

From the mathematical and physical point of view, the use of [HI] instead of [HS] brings a great benefit to our analysis. In fact its Fourier transform highlights clearly the effect of this interaction term; moreover an estimate of some properties of ε can follow.

Thus, move again to the frequency domain:

F {HI} = F {I} * F {H}
The first factor in the convolution is:

F {I} = 1 0 β(ξ) +∞ -∞ e -jωt sin Ω(ξ)t dt dξ = = -j π 2 1 0 β(ξ) [δ (ω + Ω(ξ)) + δ (ω -Ω(ξ))] dξ Let dΩ = Ω (ξ) dξ: F {I} = -j π 2 1 0 β (ξ) Ω [δ (ω + Ω) + δ (ω -Ω)] dΩ
Writing the distribution of natural frequency Ω(ξ) as the solution of the differential equation Ω = f (Ω), for an arbitrary function f , the previous integral produces the expression:

F {I} = -jω π 2 β (ω) ωf (ω)
for ω ∈ [Ω(0); Ω(1)]

F {I} = 0 elsewhere

The second factor in the convolution is:

F {H} = 1 2 δ (ω) + 1 jω
Therefore:

F {HI} = -jω π 4 β (ω) ω f (ω) + +∞ -∞ π 2 β (ζ) f (ζ) (ζ -ω) dζ
This expression provides the frequency domain counterpart of the term x * [HI]. This result shows a double contribution in the interaction term between x and the hidden variables: one is imaginary, one is real. The imaginary dissipative term is a frequency dependent damping controlled by an equivalent viscous coefficient c eq (ω) = π 4 β(ω) ωf (ω) . The function f depends on the natural frequency distribution within the hidden part of the system; more precisely the inverse of f (1/f = n/N ) is proportional to the density n of the natural frequencies over the frequency axis, so that the more the natural frequencies are dense, the more the damping effect is large.

The coefficient β is indeed:

β (Ω) = lim ω→Ω (ω -Ω) kT 0 G (ω) k0
and corresponds to the reactive part of the interaction term. For the second contribution, introduced the function Γ = πβ/2f one has:

+∞ -∞ π 2 β (ζ) f (ζ) (ζ -ω) dζ = Γ (ω) * 1 ω = π Γ (ω)
where the tilde stands, as in the previous section, for Hilbert transform. Without any further simplification, the equation of motion exhibits an interaction force that is directly the inverse transform of the previous determined term [HI]. However, accordingly with the analysis developed in the previous section, and on the same line that led to replace S by I under the hypothesis of eigenvalues "dense enough", one is pushed to estimate the determined frequency dependent terms at a given frequency ω 0 , that is the characteristic frequency around which the natural frequencies of the hidden system are located. Therefore, the imaginary part can be approximated by a viscous force represented by an equivalent viscous damper as:

c eq (ω) = π 4 β (ω) ωf (ω) ≈ π 4 β (ω 0 ) ω 0 f (ω 0 ) = c eq
The real part can be indeed roughly estimated by identifying 1/f as a delta function centered at ω 0 so that:

+∞ -∞ π 2 β (ζ) f (ζ) (ζ -ω) dζ ≈ π 2 β (ω 0 ) ω 0 = k eq
With these expressions, the time domain equation for x becomes:

m 0 ẍ + (k 00 + k eq ) x + c eq ẋ + ε (t) = f 0
For the reasons arising clearly from the previous analysis, we name again the motion of x determined by the previous equation, and for ε = 0, motion of x about ω 0 . A comparison with the result obtained in the previous section shows strong analogies, even in the form of the equivalent damping and equivalent stiffness, with the equation obtained using the Taylor expansion of G.

As a final point let us summarize the obtained result under a Hamiltonian point of view.

If a Hamiltonian system S has stiffness and mass matrices:

K = k 00 kT 0 k0 K , M = m 0 0 0 M
and lagrangian function:

L (x, ẋ) = L x (x, ẋ) + L mix (x, x) + L hid x, ẋ L x (x, ẋ) = 1 2 m 0 ẋ2 + 1 2 k 00 x 2 -f 0 x L mix (x, x) = 1 2 x kT 0 x + 1 2 X kT 0 x L hid x, ẋ = 1 2 ẋT M ẋ + 1 2 xT Kx -f T x
then, the motion of x about ω 0 (in the sense previously specified), when the other variables are unmonitored, is controlled by the new lagrangian

L x (x, ẋ) = L x (x, ẋ) + π 4 β (ω 0 ) ω 0 x 2
and by the Rayleigh dissipation function:

D ( ẋ) = π 8 β (ω 0 ) ω 0 f (ω 0 ) ẋ2
where:

β (ω 0 ) = lim ω→ω0 (ω -ω 0 ) kT 0 G (ω) k0 f (ω 0 ) = ∆ω ∆N ω0 = N n (ω 0 )
∆ω in a bandwidth around ω 0 and ∆N the number of eigenvalues of G that falls within ∆ω, or if n(ω) is the so called modal density, n(ω 0 ) is the maximum of the modal density of G. This is the fundamental result of the paper. However, as mentioned at the beginning of this section, the present approach gives the additional chance of estimating interesting properties of the remainder effect ε. The next section attacks just this problem.

Remainder term: return time and energy transfer rate

We provide in this section a nice property for ε, that is:

for t<t* the contribution of ε is negligible; t* , named return time, depends on the form of G. a. Lemma 1 As a first lemma, a simple asymptotic property for I is enlightened:

lim t→∞ I(t) = 0
This property follows from an asymptotic expansion (integration by parts) for I. Precisely (Watson's Lemma):

If Ω (ξ) = 0, ∀ξ ∈ [0, 1] :

I(t) = 1 0 β (ξ) sin Ω (ξ) tdξ = 1 t β (0) Ω (0) cos Ω (0) t - β (1)
Ω (1) cos Ω (1) t +o 1 t

If ∃ξ 0 ∈ [0, 1] : Ω (ξ 0 ) = 0:

I(t) = 1 0 β (ξ) sin Ω (ξ) t dξ = 2π t | Ω (ξ 0 )| β (ξ 0 ) cos Ω (ξ 0 ) t + o 1 √ t
This property implies that for t large enough, S(t) is well represented by the remainder r(t): S(t) ≈ r(t) , for large t that provides the behaviour of S at late time.

b. Lemma 2 Classical results in numerical integration theory, shows that:

I(t) = S(t) + r(t) r(t) = ∆ξ 2 d dξ β (ξ) sin Ω (ξ) t ξ= ξ(t)
And explicitly:

r(t) = ∆ξ 2 β (ξ) sin Ω (ξ) t + β(ξ)Ω (ξ)t cos Ω(ξ)t ξ= ξ(t) = ∆ξ 2 A(ξ) sin Ω (ξ) t + ϕ(ξ) ξ= ξ(t) A (ξ) = β 2 + β 2 Ω 2 t 2 , ϕ(ξ) = arctan βΩ β t
Thus, the late time behaviour for r may be represented by:

lim t→∞ r(t) = ∆ξ 2 β (ξ) Ω (ξ) t sin Ω (ξ) t + ϕ (ξ) ξ= ξ(t)
while it appears that:

lim t→0 r(t) = 0
This permits to conclude that the behaviour of S at early time is well represented by I:

S(t) ≈ 1 0 β (ξ) sin Ω (ξ) t dξ, at early time c.

Estimate of the return time

The results provided by the two previous lemmas can be used to give an estimate for the return time t * . The general picture for S, I and r is summarized: at early time, S and I are close because r is small, and the perturbation ε is small too; as time increases, at late time, I departs from S: in fact I vanishes and it is not for S, that is indeed close to r, and the perturbation ε becomes large. We are interested in providing the order of magnitude of the time needed for S to depart from I, that provides also an estimate of the time scale over which the perturbation ε holds small. The strategy is simple: (i) estimate the peak value S peak for both I and S reached at early time, (ii) then it is expected that I decreases (because of lemma 1), and that S initially holds close to I because r is small, (iii) at some late time t * , it is indeed expected that S increases departing from I reaching an amplitude close to its initial peak.

Let estimate S peak :

S(t) ≈ 1 0 β (ξ) sin Ω (ξ) t dξ =⇒ dS dt ≈ 1 0 β (ξ) Ω (ξ) cos Ω (ξ) t dξ ≈ 1 0 β (ξ) Ω (ξ) 1 - 1 2 Ω 2 (ξ) t 2 dξ dS dt (t peak ) = 0 =⇒ 1 0 β (ξ) Ω (ξ) dξ -t 2 peak 1 2 1 0 β (ξ) Ω 3 (ξ) dξ = 0 t peak ≈ 2 1 0 β (ξ) Ω (ξ) dξ 1 0 β (ξ) Ω 3 (ξ) dξ , S peak ≈ √ 2 1 0 β (ξ) Ω (ξ) dξ 3/2 1 0 β (ξ) Ω 3 (ξ) dξ 1/2
Let finally estimate the time t * for which S regain the value S peak ; this happens at late time, for which S can be approximated by r; evaluating all the frequency dependent functions at ω 0 :

S ≈ ∆ξ 2 β (ξ) Ω (ξ) t sin Ω (ξ) t + ϕ (ξ) ξ= ξ(t) ≈ ∆ξ 2 β (ω 0 ) f (ω 0 ) t t * ≈ 2 √ 2N f (ω 0 ) = 2 √ 2n (ω 0 )
This is also a central result of the paper. Note that the time t * has not surprisingly this expression: in fact, considering the special case in which S is a periodic function (a Fourier series), t * would be the period of S:

S = N k=1 β k sin ω k t = N k=1 β k sin 2πkt t * ω k = 2πk t * =⇒ ∆ω = 2π t * =⇒ n = 1 ∆ω =⇒ t * = 2πn
Thus, the obtained expression for t * provides essentially a proportionality between the return time and the modal density n at ω 0 . This makes clear how the interaction term in the found simplified form of an equivalent damping and stiffness is valid only when ε is small, that is up to time t * .

A final interesting result correlates the time t * to the equivalent damping and stiffness. In particular it is easy to find:

t * ≈ 4 √ 2 N k eq c eq (7)
where an intriguing proportionality is found between the return time and the equivalent damping. This result states: the faster the energy transfer from x to the hidden variables, the longer the energy storage within them. This result makes meaningful the search of special configurations of the hidden system that optimize the energy suction and storage from x: in fact if we are able to make the energy transfer as fast as possible, it comes together with the longest storage time for the moved energy.

The question is analysed in the next section.

A variational theorem for the minimum remainder term

The present section describes the conditions under which the series S(t) approaches the integral I(t) and shows that a criterion of minimum distance D(t)

D 2 (t) = C [S -I] 2 W dC = (S -I) 2 = r 2
can be satisfied with a suitable weighting function W in a prescribed domain C within a certain space Σ. A method to find such a weighting function that depends on I, which depends on ω(ξ) is described by the use of a variational approach.

The results show there exists a class of functions ω (ξ) that minimizes the distance between S(t) and I(t). In such cases, the series S(t) tends to match as close as possible the trend of the integral I(t), producing closely the apparent damping effects previously discussed, and a nearly irreversible energy transfer processes in conservative linear systems. The next section reviews the definitions and properties necessary to form the basis for the ensuing theoretical development.

The following approach is based on an average defined in a multi-dimensional space C with the use of a non trivial weighting function W. The reason for this is that, despite the difficult mathematical form, as it appears at a first glance, it can indeed lead, surprisingly, to a closed form solution of our minimization problem.

In order to correctly formulate our minimization problem, and some related constraints, let introduce some preliminary considerations.

a. Definition of average

For a set of functions s i = β i sin ω i t, i = 1, 2, . . . , N at any time t, s = [s 1 , s 2 , . . . , s N ]

T defines a point (or a vector) in the space Σ of harmonics; s exists within the hypercube

C ≡ {E × E × • • • × E}, with E ≡ [-β max , β max ], β max = max {β 1 , β 2 , . . . , β N }.
Let f (s) be an arbitrary function defined over C ⊂ Σ with the vector s ∈ Σ. In general, the average value f (t) of f (s) over C can be expressed using a weighting function P (s, I) as:

f (t) = C f (s)P (s, I) dC, dC = N k=1 ds k (8)
where the weighting function is selected to depend on s and I as described below.

As a consequence, scalar product of the two functions f (s) and g (s) in C follows as

f • g = C f (s)g(s)P (s, I) dC (9)
Similarly, the distance D(t) between f (s) and g (s) follows form:

D 2 (t) = (f -g) • (f -g) = (f -g) 2 = C [f (s) -g(s)] 2 P (s, I) dC (10)

b. Weighting Function

The weighting function P (s, I) in Eq. ( 10) is selected to have the form:

P (s, I) = N k=1 p(s k , I) ( 1 1 )
where p(s k , I) is an arbitrary function that satisfies the conditions:

Gmax -Gmax σp(σ, I) dσ = I(t) ( 1 2 ) Gmax -Gmax p(σ, I) dσ = 1 (13)
Equation ( 12) offers a comparison with the integral in Eq.( 2) for σ = β (ω) sin ωt and through a change of integration variables in Eq. ( 12) (first from dσ to dω, then to dξ): 2) and ( 14) implies that the function p(σ, I) must satisfy the following compatibility condition :

I(t) =
p(σ, I) dσ dω dω dξ = 1 (15)
The condition expressed by Equation ( 15) also implies a dependence between p (σ, I) and ω (ξ) for σ = β (ω) sin ωt.

It follows that substituting equation ( 15) in equation ( 13) yields the upper bound of ξ as : 1)can be expressed by substituting for si the averaging property -expressed for f (t)-in Eq. ( 8):

S(t) = 1 N N i=1 si = 1 N N i=1 C
s i P (s, I) dC

Further substitution for P (s, I) from Eq. ( 11) and for dC from Eq. ( 8)yields:

S(t) = 1 N N i=1 Gmax -Gmax s s i p(s i , I) ds i C (N -1) N k =i p(s k , I) ds k
The condition (13) leads the mutiplication series in the second integral to produce identity and applying condition [START_REF] Carcaterra | An entropy formulation for the analysis of energy flow between mechanical resonators[END_REF] to the first integral shows that:

S(t) = I(t) ( 1 6 )
By invoking the definition of average value in (8) leads to the fundamental result:

S (t) = C SP (s, I) dC = I(t) ( 1 7 )
provided that equation ( 15) is satisfied.

d. Constraints

The weigthing function, as it has been defined in subsection b, must satisfy the constraints ( 12) and ( 13). Equation ( 12) is automatically satisfied provided that equation ( 15) holds, as it will be used ahead. It is indeed suitable to introduce explicitly the constraint (13) in a form that makes easier to approach the following variational problem. Substituting for P (s, I) from Eq. ( 11) and for dC from Eq. ( 8) and invoking the condition in Eq. ( 13) it can be show that:

C P (s, I) dC = 1 (18)
The first derivative of Eq. ( 18) with respect to I can be expressed as: which is equivalent to stating:

∂ ∂I log P = 0 ( 19
)
that is a constraint derived from [START_REF] Carcaterra | An entropy formulation for the analysis of energy flow between mechanical resonators[END_REF] or [START_REF] Kendall | The advanced theory of statistics[END_REF]. Following the same approach, the first derivative of Eq. ( 17) with respect to I, produces a similar expression:

C S ∂ ∂I P (s, I) dC = 1 =⇒ C S 1 P (s, I) ∂ ∂I P (s, I) P (s, I) dC = 1 S ∂ ∂I log P = 1 (20)
that is a consequence of [START_REF] Nigmatullin | To the nature of irreversibility in linear systems[END_REF]. In order to find the distance between S(t) and I(t), an equivalent expression for Equation ( 20) is developed for I(t) by multiplying Eq. ( 19) by the factor I, which is independent of the integration variable s:

I ∂ ∂I log P = 0 (21)
Finally the combination of ( 20) and ( 21) produces:

(S -I) ∂ ∂I log P = 1
e. The Euler-Lagrange equation minimizing D 2 = (S -I) 2 The problem outlined at the beginning of this section 2.5 is solved indeed by a variational approach finding the minimum of the modified functional:

δ D 2 = 0 , D2 = (S -I) 2 + λ (S -I) ∂ ∂I log P
where λ is the Lagrange's multiplier to include the constraint derived in subsection d. The variation of D2 to be considered is with respect to I, that is, equivalently, with respect to ω (ξ) from which I depends, but it is easier not to consider explicitly. Therefore:

-2 (S -I) + λ - ∂ ∂I log P + (S -I) ∂ ∂I ∂ ∂I log P = 0
where the last is the Euler-Lagrange equation. This differential equation, as it can be easily verified, admits a solution of the form:

∂ ∂I log P = - (S -I) λ (22)
Condition ( 22) represents a differential equation in terms of P and its solution leads to a family of exponential functions P (s, I) = N k=1 p (s k , I). The solution to Eq. ( 22), originally obtained by Pitman and Koopman in the context of the theory of estimators, is given as:

p(σ, I) = exp {A(I)B(σ) + C(σ) + D(I)}
where A(I), B(I), C(σ), and D(σ) are arbitrary functions of their respective arguments.

Gauss function also belongs to this family of solutions and provides an excellent example that can be easily verified by substituting into equation [START_REF] Carcaterra | Pseudo-damping in undamped plates and shells[END_REF]:

p(σ, I) = 1 r √ 2π exp - 1 2 (σ -I) 2 r 2 (23)
The solution p(σ, I), with σ (β, ω), has a shape that depends on the function I(t) and on the parameter r.

Together with equation ( 23), the compatibility equation ( 15) becomes a nonlinear differential equation and its solution provides the frequency distribution ω(ξ) that minimizes the distance D:

1 r √ 2π exp - 1 2 (σ -I) 2 r 2 d dω [β (ω) sin ωt] dω dξ = 1 (24)
Equation ( 24) can be solved for ω(ξ) numerically; however, an alternative approach using density of harmonic functions, analogous to modal density in a dynamical system, produces a closed-form expression. Since dξ dω N represents the harmonic density, δ(ω), that counts the number of harmonics, or modes, contained within the frequency band dω, Eq. ( 24) directly leads to an expression for δ(ω). With dξ = 1 N dn, N being the total number of harmonics, s i , for ξ ∈ [0, 1], and dn the number of harmonics for ξ ∈ [ξ, ξ + dξ], it follows that N dω dξ ∝ dω dn = 1 δ(ω) . Substituting in equation ( 24) produces:

δ opt (t) = 1 N 1 r √ 2π exp - 1 2 (σ -I) 2 r 2 d d ω [β (ω) sin ωt] (25)
Equations ( 24) and [START_REF] Dell'isola | Piezo-ElectroMechanical (PEM) structures: passive vibration control using distributed piezoelectric transducers[END_REF] show that time appears as a parameter in the frequency distribution that minimizes the difference between S and I. Oscillators with time-dependent parameters or, equivalently, with time-varying natural frequencies, imply presence of parametrically controlled resonators or nonlinear resonators. The problem under consideration in this paper addresses linear time-invariant dynamical systems and, thus, equations ( 24) and ( 25) cannot be satisfied at all times t. Thus, the approach taken here uses the frequency distribution ω(ξ) that results from equation ( 24) or (25) for a particular time t 0 to solve the compatibility equation

p(σ, I) dσ dω dω dξ t=t0 = 1 (26)
The choice for t 0 , selection of the frequency interval [ω min , ω max ] within which ω(ξ) falls, which also depends the choice of t 0 , and the implication of their selection are discussed with examples in the next sections.

Normally, the form of Eq. ( 23) satisfies Eqs. ( 12) and ( 13) automatically for an integration domain [-∞, +∞]; however, since the actual domain is finite E ≡ [ -G max , G max ], r and I (t 0 ) must satisfy the additional constraints:

r G max , I(t 0 ) ∈ E (27)
These constraints guarantee that the function represented by equation ( 23) has its primary distribution within the interval E and therefore (approximately) satisfying equations ( 9) and ( 10).

Examples of application: set of parallel resonators

The examples given in this section illustrate application of the theory described above. Each case demonstrates how to minimize the difference between a sum of harmonic functions and the corresponding integral summation. The first example consists of a simple summation of sine functions for which β (ω) ≡ 1. In the second example, β (ω) ≡ ω represents the reaction force of a set of undamped resonators on a common rigid base. Subsections c and d examine more complex examples.

a. Simple sine series β (ω) ≡ 1 Summation of a series of N = 100 sine functions with frequencies ω i results from Eq. ( 1)by substituting for β (ω) ≡ 1:

S(t) = 1 N N i=1
sin ω i t and the corresponding integral from Eq. ( 2)becomes:

I(t) = 1 0 sin ω (ξ) t dξ
For this case, the nonlinear differential equation ( 24) becomes:

1 r √ 2π exp - 1 2 (σ -I 0 ) 2 r 2 t 0 cos ωt 0 dω dξ = 1 ( 28
)
where σ and I 0 from Eqs. ( 26) and ( 2)become

σ = sin ω (ξ) t 0 , I 0 = ξmax 0 sin ω (ξ)t 0 dξ In this case, E ≡ [-1, 1].
Restricting the selection to monotonic frequency distributions ω(ξ), so that dω dξ > 0 for ω ∈ [ω min , ω max ], implies that according to Eq. ( 28), d σ d ω = cos ω (ξ) t 0 must always be positive forω ∈ [ω min , ω max ]. It follows that assigning, for example, t 0 = π 4 , yields dω dξ > 0 for ω ∈ [0, 1]. Under these conditions, the values for I 0 may be arbitrary, except that they must satisfy the inequalities in [START_REF] Roveri | Vibration absorption using nondissipative complex attachments with impacts and parametric stiffness[END_REF] and p satisfies the conditions ( 11) and ( 15). Figures 12345illustrate the results obtained for N = 100, t 0 = π 4 , r = 0.05, ω ∈ [0, 1] and with the choice of I 0 = 0.2 and r 1, both of which satisfy equation [START_REF] Roveri | Vibration absorption using nondissipative complex attachments with impacts and parametric stiffness[END_REF]. Figure 1 represents the frequency distribution ω (ξ) determined by a numerical integration of equation ( 28) from which the set ω i (i = 1, . . . , 100) is determined by sampling 100 points equally spaced alongξ. Figure 2 represents the harmonic density and Figure 3 the time history of S(t). As shown in Figures 4 and5 for different time scales, in the time history of the series obtained using the theory developed here the strong periodicity disappears when compared with the corresponding series consisting of a linear frequency distribution with period 2πN/ω max for two different time scales.

Figures 6789show a case analogous to the previous one except for r = 0.01, representing a higher harmonic density around its peak resulting in a somewhat better performance.

The third example, shown in Figures 101112, uses t 0 = π 8 and N = 100,

r = 0.05, ω ∈ [0, 1].
The results of the examples above show that the frequency distributions satisfying the minimum distance bound requirements produce time histories that bring the summation S(t) very close to I(t), without recurrence or periodicity in its time history unlike, say, the case of a linear frequency distribution. The envelope of the summation in Eq. ( 1)decays significantly with respect to its early oscillations and without regaining its initial amplitude, following closely the same trend that its integral counterpart I(t) exhibits in Eq. (3).

b. Reaction force of a set of parallel oscillators on a rigid base, β (ω) ≡ ω Consider a set of N parallel oscillators attached to a common rigid base. Oscillators have equal mass m and uncoupled natural frequencies

ω i = k i /m,
where k i represents the stiffness of oscillator i. Impulse response of each oscillator is expressed as:

h i (t) = 1 mω i sin ω i t
The total reaction force exerted on the base by a set of N oscillators can be represented as:

S(t) = N i=1 k i h i (t) = N i=1
ω i sin ω i t S(t) has the same form as the series in (1) with G ≡ ω. In this case, for time t = t 0 , Eq. ( 26) together with Eq. ( 23) for p(σ), provides:

1 r √ 2π exp - 1 2 (σ -I 0 ) 2 r 2 [sin ωt 0 + t 0 cos ωt 0 ] d ω d ξ = 1 σ = ω (ξ) sin ω (ξ) t 0 , I 0 = ξmax 0 ω (ξ) sin ω(ξ)t 0 dξ (29)
As before, restricting the analysis only to monotonic frequency distributions ω (ξ), such that dσ dω = sin ωt 0 + t 0 cos ωt 0 > 0, and choosing, for example, t 0 = π 4 , leads to the condition that in the frequency interval ω ∈ [0, 2],

d σ dω > 0 and σ ∈ [0, 2]
. Again, r and I 0 can be assigned arbitrarily, but consistent with inequalities [START_REF] Roveri | Vibration absorption using nondissipative complex attachments with impacts and parametric stiffness[END_REF]; in this case, r = 0.1 and I 0 = 0.8. Figure 13 displays the frequency distribution obtained by solving Eq. ( 29) and Figure 14 shows the corresponding optimal modal density from equation [START_REF] Dell'isola | Piezo-ElectroMechanical (PEM) structures: passive vibration control using distributed piezoelectric transducers[END_REF]. The time history of the reaction force on the rigid base, shown in Figure 15, exhibits a rapid decay and remains at a negligibly low amplitude.

c. Pseudo-damping in conservative continuous structures, β (ω) ≡ 1 ω A continuous linear undamped dynamic system, excited by a unit impulse at point x 0 , satisfies the equation of motion:

L [w(x, t)] + m ∂ 2 w(x, t) ∂t 2 = 0
with initial conditions w(x, 0) = 0, ẇ(x, 0) = δ(xx 0 ), where δ is the Dirac's distribution. L [ ], w(x, t), m represent the system operator, the displacement response and the mass density, respectively. The general response of such a linear system can be represented by its orthonormal modes Φ i (x) and principal co-ordinates q i (t), as

w(x, t) = N i=1 Φ i (x), q i (t) q i (t) = A i sin ω i t , A i = m ω i Φ i (x 0 )
Then its impulse response at x 0 can be represented by the series expression S(t) in Eq. ( 1), with G i = m ωi Φ 2 i (x 0 ). In general, in the absence of damping, this finite series, a superposition of pure sine functions, exhibits an almost-periodic trend. For example, the case of a Fourier series of sine functions with linearly distributed frequencies ω i = iω 0 , where ω 0 is the fundamental frequency, becomes periodic. As before, a decaying trend in S(t) is expected only in the presence of energy dissipation. However, as shown in previous studies that in cases where condensation points exist within the frequency distribution or, equivalently, natural frequencies accumulate around a particular frequency, impulse response of that linear system exhibits a decaying characteristic even in the absence of dissipation sources, a phenomenon referred here as near-irreversibility or apparent-damping.

Application of the theory developed in this paper to the continuous system described above provides a theoretical basis to the numerically obtained results in earlier studies and demonstrates how a class of frequency distributions ω i can produce apparent-damping.

As an example, consider a simply-supported beam as a prototypical lin- 1), and retaining only the odd terms:

ear system for which Φ i (x 0 ) = 2 m L sin πix0 L , x0 L = 1 2 . Substituting for G i = 1 ωi 2 L sin πi 2 2 in Eq. (
S(t) = N/2 i=1 2 L 1 ω 2i-1 sin ω 2i-1 t
Selecting t 0 = π 4 , the function σ = 1 ω sin ωt 0 has a monotonically increasing trend, for example, within the interval ω ∈ [START_REF] Carcaterra | Transient energy exchange between a primary structure and a set of oscillators: return time and apparent damping[END_REF][START_REF] Murat Koç | Energy sinks: vibration absorption by an optimal set of undamped oscillators[END_REF]. Choosing values I 0 = 0.05, r = 0.005, which satisfy the conditions [START_REF] Roveri | Vibration absorption using nondissipative complex attachments with impacts and parametric stiffness[END_REF], the frequency distribution and the corresponding modal density can be obtained from equations [START_REF] Dell'isola | Passive damping of beam vibrations through distributed electric networks and piezoelectric transducers: prototype design and experimental validation[END_REF] and [START_REF] Dell'isola | Piezo-ElectroMechanical (PEM) structures: passive vibration control using distributed piezoelectric transducers[END_REF], as shown in Figures 16 and17. The impulse response of the beam with such a frequency distribution is shown in Fig. 18 for N = 200 (but includes only the 100 odd modes). The impulse response shows a rapid decay reminiscent of the impulse response of a highly-damped system, being indeed the system is Hamiltonian. d. Pseudo-irreversible energy transfer between a single dof resonator and a set of parallel oscillators, β (ω) ≡ ω 3 Figure 19 depicts the system under consideration in this section, which consists of set of resonators with natural frequencies ω i (i = 1, 2, . . . , N), that are connected in parallel to a common principal structure. The system does not possess any means of energy dissipation. For a very large number of oscillators N, approaching infinity, but with a constant total mass, the attached oscillators can be considered as a continuous distribution with the equations of motion for the coupled system expressed as:

       M ẍM (t) + K M x M (t) + 1 0 k(ξ) (x M (t) -x(ξ, t)) dξ = 0 mẍ(ξ, t) -k(ξ) (x M (t) -x(ξ, t)) = 0 (30)
where M , K M , x M are the mass, stiffness and displacement of the master structure, respectively; m, k, x represent the same quantities of the distributed oscillators in the attached set.

Several studies have shown that such a distribution of oscillators, produce a damping effect on the principal mass (1; 2; 3; 4; 5; 6; 7; 8) as N approaches infinity.

An alternative derivation of this result, presented in the Appendix, shows that the impulse response of the principal oscillator progressively decays and vanishes asymptotically with time. Energy initially imparted to the principal structure migrates to the attached set of infinite number of oscillators that have frequencies that fall within a finite bandwidth, where it remains indefinitely. As discussed earlier, it is commonly accepted that, in general, such irreversible energy transfer does not hold for a finite N (8). However, as the following application of the theory developed in this paper shows, there exist particular frequency distributions which afford a nearly-irreversible energy transfer even for a finite set of oscillators.

Considering the second of equations ( 29), the displacement of the continuous set of resonators in terms of the master response can be expressed by the convolution integral:

x(ξ, t) = ω n (ξ) t 0 x M (τ )H (t -τ ) sin ω n (ξ) (t -τ ) dτ
where H is the Heaviside function. Introducing this expression into the first of Eq. ( 29), an integro-differential equation results in terms of x M :

M ẍM (t) + K M x M (t) + x M (t) 1 0 k(ξ) dξ+ - t 0 x M (τ ) 1 0 mω 3 n (ξ)H(t -τ ) sin ω n (ξ)(t -τ ) dξ dτ = 0
which can be also expressed as:

M ẍM (t) + (K M + k)x M (t) -x M (t) * [I(t)H(t)] = 0 (31)
where I(t)H(t) is the kernel of the integral part of the previous equation and

k = 1 0 k (ξ) dξ, I(t) = 1 0 mω 3 n (ξ) sin ω n (ξ) t dξ (32)
In the case of a finite set of N resonators, the equation of motion takes a different form where integrals over ξ are substituted by summations. Thus, Eq. ( 31) remains applicable provided that k =

N i=1 k i and I(t) is replaced by its discrete counterpart S(t) = 1 N N i=1 mω 3 i sin ω i t: M ẍM (t) + (K M + k)x M (t) -x M (t) * [S(t)H(t)] = 0 (33)
The apparent damping and near irreversibility as manifested by the decay characteristics of the impulse response result from the application of the present theory by considering G (ω) = mω 3 with σ = m ω 3 sin ωt 0 and searching for the optimum frequency distribution.

As an example, consider a master structure, with an uncoupled natural frequency ω M = 1, with N = 100 attached oscillators. Assuming t 0 = π 4 and searching for a monotonic frequency distribution ω (ξ), it follows that, within the frequency intervalω ∈ [0, 2], dσ dω > 0 and σ ∈ [0, 8]. The values of r and I 0 (r = 0.4 and I 0 = 0.6) are selected to be consistent with inequalities [START_REF] Roveri | Vibration absorption using nondissipative complex attachments with impacts and parametric stiffness[END_REF], and to assure that the function represented by equation ( 23) has its peak around ω = ω M = 1.

Figures 20 and21 show the frequency distribution and the frequency density of the attached oscillators determined by solving equations ( 24) and [START_REF] Magionesi | Insights into the energy equipartition principle in undamped engineering structures[END_REF]. Figure 22 shows the master response following an impulse applied at t = 0, which illustrates how a significant part of its energy is transferred to the set of oscillators and remains there without returning back to the master, producing an irreversible energy transfer.

Frequency (rad/s) ( ) Frequency (rad/s)

!(») » Figure 6
Harmonic density (s/rad) ( ) Frequency (rad/s)

S(t)

!(») » Figure 10
Frequency (rad/s)

!(») » Figure 13
Harmonic density (s/rad) ( )

S(t) t

Frequency (rad/s) Frequency (rad/s)

!(»)

!(»)

UNISAT: An Engineering Application

This section presents an application of a vibration damper based on a cluster of beams, based on the theory developed in section B. The theory is applied to a cluster of continuous beams structures attached to a continuous master, outlining an analysis for predicting the expected performances of the damper with particular emphasis to an optimal design of the device.

A built-up system is applied to a satellite aerospace structure to be launched next year, capable of absorbing the vibration energy at the liftoff, accordingly with the findings of the theory outlined in the paper. The experimental results illustrate the feasibility and the attractiveness of this new damping technique. a. Pseudo-damping effect induced by a cluster of beams Aim of the section is to derive the master-cluster coupled equations of motion. The cluster consists of a set of N beams, all clamped at the same support AB, attached to the main structure (master) by the connector QP as illustrated in figure 1. AB and QP are assumed to be rigid (their stiffness much higher than those of the beams) and the attachment between QP and the master is assumed to be a point connection.

The coupled equations of motion The master, a continuous linear undamped structure, is described by the governing equation:

L [w(x, t)] + m M ∂ 2 w(x, t) ∂t 2 = Rδ(x -x P ) ( 3 4 )
where δ is the Dirac's distribution, L [ ], w(x, t), m M are the structural operator, the structural displacement and the mass density. x is the space variable along the master and R is the reaction force exerted by the cluster on the master at x P .

The master is coupled to the beams of the cluster through R:

R = N i=1 T i = N i=1 B i ∂ 3 w i ∂ξ 3 i ξi=0 (35)
where T i , B i , w i , ξ i are: the shear force at the clamping support, the bending stiffness, the flexural displacement (in the same direction of w) and the abscissa along the i-th beam of the cluster, respectively. The N + 1 coupled equations of motion for the N beams and for the master are:

           L [w(x, t)] + m M ∂ 2 w(x, t) ∂t 2 -δ(x -x P ) N i=1 B i ∂ 3 w i ∂ξ 3 i ξi=0 = 0 B i ∂ 4 w i (ξ i , t) ∂ξ 4 i + m i ∂ 2 w i (ξ i , t) ∂t 2 = m i ∂ 2 w(x P , t) ∂ t 2 , i= 1, 2, . . . N (36)
where m i is the mass per unit length of the i-th beam, and the beams of the cluster are forced by the inertial term related to the motion of AB that moves as P . Φ ih (ξ i ) and q ih (t) are the i-th orthonormal mode of the h-th beam (clamped-free) and the associated principal co-ordinates of the master, respectively. The beam vibration fields within the cluster are described by:

w i (ξ i , t) = ∞ i=1 Φ i h (ξ i )q i h (t) ( 3 7 )
Substitution for these expressions into the beams equations (the second of system [START_REF] Cercignani | Ludwig Boltzmann-The Man Who Trusted Atoms[END_REF]) and the use of the orthonormality conditions produces:

qih (t) + ω 2 ih q ih (t) = L ih ∂ 2 w (x P , t) ∂t 2 (38)
where ω ih is the i-th natural frequency of the h-th beam and

L i h = li 0 m i Φ i h (ξ i ) dξ i (39)
For equations [START_REF] Dell'isola | Continuum modelling of piezoelectromechanical truss beams: An application to vibration damping[END_REF], the convolution-form solutions hold:

q ih (t) = Γ(t) L i h ω i h sin ω ih t * ∂ 2 w (x P , t) ∂t 2 (40)
where Γ(t) is the Heaviside distribution. Equations ( 37) and ( 40) permit, after simple mathematics, to rewrite the equation of the master (first of equations ( 36)) in the form:

L [w(x, t)] + m M ∂ 2 w(x, t) ∂t 2 -δ(x -x P ) [S(t)Γ(t)] * ∂ 2 w (x P , t) ∂t 2 = 0 (41)
where

S(t) = ∞ h=1 N i=1 U ih ω i h sin ω ih t, U ih = B i L ih ∂ 3 Φ i h ∂ξ 3 i ξi=0 (42)
Equation ( 41) is integral-differential in terms of the master displacement w. This prototype equation has been studied in the previous section and remarkable properties of the kernelS(t)H(t) have been established. Of particular interest, the possibility, under certain conditions, of replacing the summation through the set of beams (index i) appearing into S(t) by an integral, enlightening in this way the damping effect the cluster produces in the master motion. Physically, this means that the beams of the cluster, having a discrete spectrum of natural frequencies within the set, are replaced by a continuous distribution of frequencies with a spacing tending to zero. Introduction of the dummy variable χ, varying between 0 and 1 through the set of resonators with ∆χ = 1/N , allows summation [START_REF] Dell'isola | Almansi-type boundary conditions for electric potential inducing flexure in linear piezoelectric beams[END_REF] to be written as an integral I(t):

S(t) = ∞ h=1 N N i=1 U h (χ i ) ω h (χ i ) sin ω h (χ i )t 1 N ≈ ≈ I(t) = ∞ h=1 N 1 0 U h (χ) ω h (χ) sin ω h (χ)t dχ (43)
whereω h (χ) represents the continuous distribution of the h-th natural frequency through the set of the beams that depends on the section and lengths of the beams within the cluster. This integral can be fruitfully manipulated replacing the integration variable χ byω h . Let dω h (χ) dχ = f (ω h ) and change the integration limits:

I(t) = ∞ h=1 N ω h (1) ω h (0) U h (χ) ω h f (ω h ) sin ω h t dω h (44)
f (ω h ) is the frequency density of the h-th mode of the beams. The final step consists in expressing U h (χ) in terms of ω h . As shown in Appendix A, for a set of beams of different lengths and cross section area:

U h (χ) = 4B (χ) β 2 h l 3 (χ) sinh β h -sin β h cos β h -cosh β h 2 (45)
where β 1 = 1.875, β 2 = 4.694, β 3 = 7.855 . . . etc. The distribution of natural frequencies ω h (χ) depends on the bending stiffness B(χ) and the length l(χ) distributions through ω h (χ) =

β 2 h l 2 (χ) B(χ) m (χ)
, where m (χ) is the mass per unit length across the beams set. Substitution for this last expression into equation [START_REF] Alessandroni | Optimal piezo-electromechanical coupling to control plate vibrations[END_REF] gives:

U h (χ) = µ h ω 2 h (χ) m (χ) , µ h = 4 β 2 h sinh β h -sin β h cos β h -cosh β h 2
where m(χ) = m (χ)l (χ) is the beams mass distribution within the set; m(χ) can be replaced by m(ω h ), since, through the distribution ω h (χ), a frequencyω h remains associated at anyχ. Thus, S can be conveniently approximated by I(t):

S(t) ≈ I(t) = ∞ h=1 N ω h (1) ω h (0) µ h ω h m (ω h ) f (ω h ) sin ω h t dω h (46)
and the final equation for the master motion is obtained by substitution of ( 46) into [START_REF] Maurini | Comparison of piezoelectronic networks acting as distributed vibration absorbers[END_REF]:

L [w(x, t)] + m M ∂ 2 w(x, t) ∂t 2 + -δ(x -x P )    ∞ h=1 N Γ(t) ω h (1) ω h (0) µ h ω h m(ω h ) f (ω h ) sin ω h t dω h    * ∂ 2 w (x P , t) ∂t 2 = 0
(47) The third addend represents the global action of the cluster on the master.

Equivalent damping of the cluster The crucial point of the present theory relies in showing how the frequency representation of (47), through Fourier transform F ,

L [W (x, Ω)] -m M Ω 2 W (x, Ω)+ -δ(x -x P )F         ∞ h=1 N Γ(t) ω h (1) ω h (0) µ h ω h m(ω h ) f (ω h ) sin ω h t dω h    * ∂ 2 w (x P , t) ∂t 2      = 0
(48) produces for the third term a nonzero imaginary part, implying the cluster effect amounts to an equivalent dissipation.

The third addend can be indeed written as:

W (x P , Ω) ∞ h=1 -Ω 2 N F {Γ(t)} * F      ω h (1) ω h (0) µ h ω h m(ω h ) f (ω h ) sin ω h t dω h      (49)
where:

F {Γ} = 1 2 δ (Ω) + 1 jΩ (50)
and

F      ω h (1) ω h (0) µ h ω h m(ω h ) f (ω h ) sin ω h t dω h      = = ω h (1) ω h (0) µ h ω h m(ω h ) f (ω h ) +∞ -∞ e -jΩt sin ω h t dt dω h = = -jπ ω h (1) ω h (0) µ h ω h m(ω h ) f (ω h ) [δ (Ω + ω h ) -δ (Ω -ω h )] dω h = = -jπ µ h Ωm(Ω) f (Ω) Π h (Ω) (51)
with

Π h =    1 for Ω ∈ [ω h (0), ω h (1)
] 0 elsewhere 50 Therefore:

F {Γ(t)} * F      ω h (1) ω h (0) µ h ω h m(ω h ) f (ω h ) sin ω h t dω h      = = -jΩ π 2 µ h m(Ω) f (Ω) Π h (Ω) - +∞ -∞ πµ h Π h (ζ) ζm(ζ) f (ζ) (ζ -Ω) dζ
In the integral on the rhs, the kernel of the Hilbert transform H is recognized:

+∞ -∞ πµ h Π h (ζ) ζ m(ζ) f (ζ) (ζ -Ω) dζ = π 2 µ h H Π h (Ω) Ωm(Ω) f (Ω)
Thus, equation (48) takes the final equivalent form:

L [(x, Ω)] -m M Ω 2 W (x, Ω)+ + δ(x -x P )W (x, Ω)N ∞ h=1 -jΩ π 2 µ h Ω 2 m(Ω) f (Ω) Π h (Ω) + -π 2 µ h Ω 2 H Π h (Ω) Ωm(Ω) f (Ω) = 0
Neglecting the real part of the series (small compared to L [W (x, Ω)]m M Ω 2 W (x, Ω)), the equation of the master finally reads:

L [W (x, Ω)] -m M Ω 2 W (x, Ω)+ + jΩ π 2 N Ω 2 m(Ω) f (Ω) ∞ h=1 µ h Π h (Ω) W (x, Ω)δ(x -x P ) = 0
This expression provides the frequency domain counterpart of the master equation and apparently suggests that the set of oscillators introduces a frequency dependent equivalent damping

C eq (Ω) = π 2 N Ω 2 m(Ω) f (Ω) ∞ h=1 µ h Π h (Ω) (52)
even in the absence of any dissipation effect.

An alternative and still significant expression of the apparent damping is obtained. In fact dm = m(Ω)N dχ = m(Ω) dN is the mass of the number dN of beams having natural frequency within the bandwidth [Ω, Ω + dΩ]. Thus, equation (52) becomes:

C eq (Ω) = π 2 Ω 2 dm dΩ ∞ h=1 µ h Π h (Ω) (53)
b. Physical considerations and properties of C eq Physical meaning of C eq : frequency and time domain The presence of a damping term in the master equation, in the absence of any energy dissipation within the whole system, amounts to an effect of energy transfer from the master to the beams, equivalent for the master to a loss of energy.

In this view the cluster works as a vibration or shock absorber. However, deriving C eq , in equation ( 43) an approximation is made replacing the summation S(t), through the set of beams, by the integral I(t), and the simple frequency domain expressions (52) or (53) are produced. These hold only if suitable conditions for this replacement occur, and the question can be more conveniently examined in time domain. In general S and I(t) differ indeed for a remainder term R = I -S. As shown in (3), R is small at early times and, within a suitable time interval [0, t * ], R ≈ 0. Thus, for t smaller than t * , S ≈ I, |R| |S| ≈ |I|, the replacement of S by I is legitimated and equations ( 52) and (53) hold.

More in detail, applying the argument based on the wave analogy illustrated in (3), for the present case, t * can be estimated by (see Appendix B):

t * ≈ 2πN f (Ω tun ) (54)
where Ω tun is the frequency in correspondence of which the set of beams is tuned. It follows that for t ∈ [0, t * ], C eq is the actual apparent damping of the master, while for later times the cluster effect is not anymore correctly represented by it. Physically, after t * , the energy temporarily released to the set of resonators is indeed suddenly returned back to the master. In section 2 it is shown that, in some remarkable cases, this energy return phenomenon can be prevented. In fact, R depends on the functions ω h (χ) and a special class of them, able to minimize a suitable average of R, is found, permitting the replacement of the summation by an integral with a minimized error. In a sense, these 'optimal' frequency distributions makes t * the largest possible, and the energy return to the master is not practically observed. In fact, a large t * , in real structures, is substantially equivalent to an infinite t * . Indeed even small dissipation transforms the most of vibration energy trapped in the cluster into heat and the energy returned back to the master, after t * , is only a negligible fraction of that trapped initially.

Finally, as shown in the previous section 2, the return time and apparent damping are related: the faster the energy transfer from the master to the cluster, i.e. the higher C eq , the longer t * , as it also appears by comparison of equations ( 53) and (54), from which results C eq (Ω M ) ∝ t * .

Thus, as detailed in the next section, an effective design of the damper gains an important benefit in using these optimal frequency distributions, having better performances both in terms of C eq and t * . c. Optimal frequency distributions and properties of C eq (Ω) Theoretically, the best frequency distribution ω h (χ) is found in found in the previous section by solving for a functional minimization of the remainder square average. The found solutions belong to exponential families. Among them, that proposed ahead is rather simple and effective:

ω opt (χ) = W 2 2 χ -1 |2χ -1| e α|2χ-1| -1 e α -1 + 1 + ω min
where χ ∈ [0, 1], W = ω maxω min and α is a parameter that modifies the shape of the frequency distribution (its optimal value for best energy absorption is 2.5). Assuming ω 1 (χ) = ω opt (χ), i.e. the optimal distribution across the set for the first mode, the set of the beams lengths within the cluster follows as l

(χ) = β1 √ ωopt(χ) 4 B
m , assumed that all the beams of the cluster have the same bending stiffness B and same cross section area (thickness h, width b, mass per unit length m ) differing only for their lengths.

Thus, the frequency distributions for the h-th mode across the beams set:

ω h (χ) = β h β 1 2 ω opt (χ) ( 5 5 )
All these distributions present a typical trend forχ ∈ [0, 1]: a rapid growth close toχ = 0, an almost flat branch aroundχ tun = 0.5, and again a rapid growth close to χ = 1, as shown in figure 2. Consequently, 1 f = 1 dω h /dχ is a Gaussian-like function with sharp peak aroundχ tun , to which corresponds the tuning frequency the values around which the natural frequencies of the distribution ω h (χ) condense (see figure 3). An intuitive understanding of the optimal characteristics in terms of t * and C eq of these distributions is seen from expressions (52) and (54) that contain the factor 1 f , implying they confer a large t * and a large C eq around the tuning frequency ω tun h .

ω tun h = β h β 1 2 ω max + ω min 2 ,
Besides these typical peaks of the equivalent damping related to 1 f , expressions ( 52) and (53) reveal other remarkable properties of C eq discussed ahead.

Figure 3. Gauss-like frequency density distribution within the cluster.

From equation (54) the frequency average damping can be estimate. The frequency intervals [ω h (0), ω h (1)] can be partially overlapped or not. In the last case (as in the built-up device here presented) ω h (1) < ω h+1 (0) for any h, meaning β h+1 β h lmin lmax > 1, where l min , l max are the minimum and maximum beams lengths through the set. Thus, in each frequency window Π h , the damping expression simplifies as C eq h (Ω) = π 2 Ω 2 dm dΩ µ h . An estimate of its average over the frequency bandwidth Ω ∈ [ω h (0), ω h (1)] is:

Ceq h = 1 W h ω h (1) ω h (0) π 2 Ω 2 dm dΩ µ h dΩ ≈ πω 2 tun h 2W h µ h M cl (56)
where W h = ω h (1)ω h (0) and M cl are the h-th device tuning bandwidth and the total mass of the cluster, respectively. Finally a relationship between M cl and N holds:

M cl = m N i=1 l (χ i ) = m N N i=1 l (χ i ) ∆χ ≈ m N 1 0 l (χ) dχ = = m N β 1 4 B m 1 0 1 ω opt (χ) dχ (57)
Equations ( 52), ( 54), ( 55), ( 56) and (57) provide the basis for the cluster design.

d. Cluster design and performances of the built-up device A device based on the previous theory has been designed to be used on board of UNISAT, that stands for UNIversity SATellite, a permanent space project developed at the University of Rome La Sapienza by the Gauss Group. It is a small scientific satellite (14 kg-20 kg depending on the flying payload) launched four times in orbit (2000)(2001)(2002)(2003)(2004)(2005)(2006) and next launch, equipped with the presented vibration suppressor, is planned in 2008. Severe vibrations occur to the electronic equipment of the satellite during the lift-off operations of the carrier and the present device is aimed at reducing shock and vibration on the plate carrying the electronic package. The material used for the damper is still (namely ρ = 7780 kg/m 3 , E = 187.5 GPa,) and the cluster of beams are realized by milling machines from a still sheet of thickness h = 0.6mm. The maximum allowed room on board for the device is 90 mm × 90 mm × 40 mm, with a maximum allowed mass agreed by the satellite designers equal to 150 g.

The design procedure follows the steps ahead:

• a) Frequency bandwidth and tuning frequency. • c) Lengths distribution. The length distribution (see figure 6) is determined on the basis of ω opt (χ) as

l (χ) = β1 √ h √ ωopt(χ) 4 E
12ρ , that yields minimum and maximum lengths of 30.7 and 35.6 mm, respectively. Considering that the device is made of symmetric beams (see figure 6), the maximum dimension along the beams axes is twice the maximum beam length, i.e. 71.2 mm, to which the clamped part of width of the beams, 10 mm, should be added leading to a maximum dimension equal to 81.2 mm that satisfies the design constraint (max 90 mm).

• d) Number of beams. The number N of the beams is to guarantee a return time large enough to prevent energy comes back to the master. As it appears indeed form equation (56), once given the frequency bandwidth and the tuning frequency, only the total mass of the cluster has an effect on C eq , but not the number of beams over which the mass is spread.

t * is the time the energy is stored within the beams. It is desired that only a small fraction r of the initially trapped energy should be returned to the master. If η is the hysteretic damping factor of a beam, roughly energy decays as e -ηωt for a harmonic motion at frequency ω. Lower frequencies have a slower decay, and the lowest decay is with the lowest frequency ω min within the cluster. Thus, we can ask the factor e -ηωmint * be equal to r: e -ηωmint * = r =⇒ t * = -ln r ηω min and from (54): N = -f (ω tun ) ln r ηω min (58) that provides the order of magnitude for the number of beams to be included within the cluster. Since f (ω tun ) ≈ 500 rad/s, ω min ≈ 2350 rad/s, η ≈ 0.03 and assuming r = 0.01, N ≈ 30.

• e) Width and gap.

The maximum allowed width of the device is D = 90 mm, g is the gap between two adjacent beams, and b the width, then:

D = N b + (N -1)g (59)
Since, as it appears from equation (56), it is convenient to increase the mass of the cluster, i.e. b must be the largest, compatibly with equation (59). This implies g must be the smallest. Actually, it depends on the smallest thickness of the milling cutter disk that generates the gap between the beams, in this case is g = 1 mm. Thus, form the previous equation follows b = 2 mm.

• f) Maximization of the cluster mass (multiple layer).

The previous analysis completely defines the optimal cluster satisfying the prescribed design requirements. However, since the maximum allowed height of the whole device is up to 40 mm, a multiple layered structure can be hosted, with the advantage of increasing the total mass of the cluster, as suggested by equation (56). Thus, three equal sets of beams are in column, separated by two aluminum spacers with height 5mm that guarantee each beam does not hit the upper or the lower beam when undergoing the maximum accelerations at lift-off provided by the satellite builder.

The final mass M DEV of the complete assembled device (represented in figure 7) is 130 g. Finally, the validation experiments are performed following the lines explained below.

As a first step, the best location for attachment point P is identified (fig. 9). An electro-dynamic shaker excites the structure with a spectrum similar to that meet in operating conditions and the maximum displacement point is identified as P. At the same time the drive point FRF at P is determined and a peak frequency identified (see figure 4, frequency peak at 440 Hz) used for determining bandwidth and tuning frequency.

A first comparison is made between the theoretical C eq (Ω) and the experimentally identified apparent damping of the device. The FRF at the device attachment point is determined and the experimental apparent damping C EXP eq (Ω) is identified by fitting the obtained response by using a theoretical FRF of the form:

F RF T H (Ω) = 1 -Ω 2 M DEV + jΩC eq (Ω)
where M DEV is known. The plot of C EXP eq (Ω) is shown in figure 10, and the comparison with figure 8 shows a good agreement. The smoothness of the theoretical curve is because the developed theory uses an integral instead of the discrete summation through the set, while the experimental peaks corresponds to the separate resonance frequencies of the single beams of the cluster.

The built-up device is then installed on the satellite plate. Comparison of the new attenuated drive point FRF and the old one is shown in figure 12 Frequency (Hz) Apparent damping coefficient (kg/s) in the frequency band covered by the first mode of the cluster of beams. In figure 13, the second modes of the device are indeed within the frequency bandwidth 2200-3000 Hz. It appears how also in this higher frequency band the damper effectively reduces the amplitude of vibration.

Frequency (Hz)

  the interval [ω min , ω max ], σ = β (ω) sin ωt is single-valued and σ ∈ [-β max , β max ]. A comparison of equations (

  max yielding ξ max = 1, which leads to the conclusion about the bounds of ξ as: ξ ∈ [0, 1]. c. Average of S The average of the function S(t) in Eq. (

  P (s, I)] P (s, I) dC = 0

Figure 1

 1 Figure 1

Figure 2 35S

 2 Figure 2

Figure 4 36S

 4 Figure 4

Figure 5

 5 Figure 5

Figure 8

 8 Figure 8

  t

Figure 9

 9 Figure 9

Figure 16 Harmonic

 16 Figure 16

Figure 19

 19 Figure 19

Figure 20 Harmonic

 20 Figure 20

Figure 1 .

 1 Figure 1. Sketch of the mastercluster coupling.

Figure 2 .

 2 Figure 2. Typical optimal frequency distribution within the cluster.

Figure 6 .

 6 Figure 6. Optimal length distribution within the set of beams of the device.

Figure 7 .

 7 Figure 7. View of the build up device tuning frequency 440 Hz, total weight 130 g.

Figure 8 .

 8 Figure 8. Apparent damping coefficient vs frequency of the device theoretically determined with an infinite number of resonators

Figure 9 .

 9 Figure 9. Selection of the test "point P " on board at which the suppression vibration device is applied.

Figure 10 .

 10 Figure 10. Experimental measurement of the apparent damping of the device (tuning frequency 438 Hz).

Figure 11 .

 11 Figure 11. View of the final installation of the vibration suppression device on board of UNISAT.

Figure 12 .

 12 Figure 12. Comparison of the FRFs at "point P" with (black curve) and without (circled curve) the vibration suppression device.

Harmonic density (s/rad) ( ) Frequency (Hz) Spectrum of the response ATP On the basis of the experimental drive point frequency response at the attachment point P for Ω ∈ [0, 1000] Hz (see figure 4), the frequency bandwidth W and its central frequency ω tun are selected. Looking at the highest peak, the choice is ω tun = 440 Hz, W = 140 Hz.

• b) Optimal frequency distribution.

With ω tun = 440 Hz, W = 140 Hz, ω max , ω min are determined, and they completely define the optimal frequency distributions given by equation (55) (see figure 5) and the related values ω tun h , W h .