Open access

Rotordynamic Stabilization of Rotors on Electrodynamic Bearings

Written By

J. G. Detoni, F. Impinna, N. Amati and A. Tonoli

Submitted: March 8th, 2012 Published: October 2nd, 2012

DOI: 10.5772/51350

Chapter metrics overview

1,945 Chapter Downloads

View Full Metrics

1. Introduction

In the last decades the advance in the semiconductors technology for power electronics has dictated a growing interest for high rotational speed machines. The use of high rotational speeds allows increasing the power density of the machine, but introduces some critical aspects from the mechanical point of view. One of the most critical issues to be dealt with is the difficulty in operating common mechanical bearings in this condition. For this reason alternatives for classical ball and roller bearings must be found. In this context, active magnetic bearings represent an advantageous alternative because they are capable of supporting the rotating shaft in absence of contact. Nevertheless, the high cost associated with this kind of system reduces their applicability.

A promising system for supporting high rotational speed machines in absence of contact and with relatively low costs, widening the range of applications, is the electrodynamic suspension of rotors [1], [2], [3], [4], [5]. Systems capable of realizing this concept are commonly referred to as electrodynamic bearings (EDB). They exploit repulsive forces due to eddy currents arising between conductors in motion relative to a magnetic field. The supporting forces are generated in a completely passive process, thus representing an increase in the overall reliability of the suspension with respect to active magnetic bearings. Nevertheless, electrodynamic bearings have drawbacks. The eddy current forces that provide levitation produce an energy dissipation that may cause negative damping resulting in rotordynamic instability.

Because the rotor may present an unstable behavior, it is necessary to study the dynamic response of the suspension in order to guarantee stable operation in the working range of speed. This can be achieved by introducing nonrotating damping in the system, but the choice of the damping elements is not obvious, requiring an accurate modeling phase. The present paper presents the development of a dynamic model of the entire suspension that is used to study the mechanical properties of the supports that allow guaranteeing rotordynamic stability. A simple optimization procedure is used in order to identify the characteristics of an elastic support placed in between the electrodynamic bearing’s stator and the casing of the machine. The use of anisotropic supports to improve the stabilization characteristics is also investigated, and optimal conditions are identified.


2. Dynamic model of EDBs

To describe the dynamics of the eddy currents inside the coils and also the dynamic effects of electrodynamic bearings on rotors supported by them, we make the assumption that the rotor rotates at constant angular speed Ω ( θ = Ω t ). Assuming constant or slowly varying rotational speed is commonly done in rotordynamics[6] and does not reduce the validity of the model.

The systems under analysis are shown in Figure 1a and Figure 1b. The first presents a schematic representation of a heteropolar EDB with the magnetic field generated using a two pole pair Halbach array. The second is a scheme of a homopolar EDB having a radial flux configuration. Different configurations are possible and can be studied using the same models presented in this paper.

Figure 1.

Scheme of possible configurations of electrodynamic bearings.(a) Heteropolar configuration; (b) homopolar configuration.

To write an equation that describes the behavior of the current in the electric circuit of the coils the electric circuit where the current flows must be defined. Figure 2a presents the electric circuit where the terms R c and L c are the resistance and inductance of the coil. For some applications it may be interesting to connect inductive loaded circuit in series with the coil [1], [2]. For this reason the terms of a generic passive shunt R a d d and L a d d are introduced in the model. The mutual inductance between the coils has been neglected. The orthogonality between the coils justifies this assumption.

Figure 2.

Schemes of the EDB’s rotor describing the main model’s variables. (a) Geometrical variables; (b) Electric circuit of the rotor’s coils taking into account of a generic passive shunt.

2.1. Eddy currents and bearing’s forces

Developing the modeling of the electrical equations in terms of complex quantities allows a strong simplification of the system’s equations. Hence we define the main geometrical and electrical variables described in Figure 2a and Figure 2b and necessary for the modeling. The Lagrangian coordinate representing the displacement of the geometric center of the rotor C relative to the axis of the magnetic field O in the Cartesian reference frame is given by q c in the form:

q c = x + j y E1

The electric current inside the coils in complex coordinates is written as:

i = i 1 + j i 2 E2

Considering the above defined variables, the state equation describing the dynamics of the eddy currents inside the coils can be written as:

d i d t = λ ˙ L R L i E3

where λ is the magnetic flux generated by the permanent magnets and linking the rotor’s coils. The expression describing the flux linkage can be expressed in complex coordinates as:

λ = q c Λ 0 e j( p - 1 ) Ω t E4

In this expression the term Λ 0 is a coefficient that gives the variation of magnetic flux linkage due to a lateral displacement of the rotor and has units of Wb/m whereas p represents the number of pole pairs of the magnetic field and is given by an even number ( p = 0 , 2 , 4 , ). Notice that a magnetic field where p = 0 is a homopolar magnetic field while p 0 gives place to heteropolar magnetic fields.

Combining Eq. (3) and Eq. (4) allows writing the system of equations coupling the rotor’s motion and the induced current as:

d i d t = Λ 0 L ( q ˙ c + j ( p 1 ) q c Ω ) e j( p 1 ) Ω t R L i F q = i Λ 0 e -j( p 1 ) Ω t E5

In the equation the current i is the state variable and the bearing’s reaction force F q is the output equation. Although linear, this equation has periodically time-varying coefficients. Performing a change of variables, substituting the state variable i with the output variable F q it is possible to obtain a set of equations having constant coefficients as:

F ˙ q = Λ 0 2 L ( q ˙ c + j ( p 1 ) q c Ω ) F q ( R L + j ( p 1 ) Ω ) i = F q Λ 0 e j ( p 1 ) Ω t E6

This set of equations allows calculating the reaction forces generated by EDBs of both homopolar and heteropolar configurations, and can be used to study the dynamics of rotors supported by EDBs.


3. Jeffcott rotor on EDBs

Due to the nature of the phenomena, studying the dynamics of a rotor on magnetic bearings requires one to consider that the center of the rotor is moving relative to the stator. In the specific case of the electrodynamic bearing, this means that the center of the conductor (point C) is moving relative to the magnetic field (point O). Equation (6) takes this into account. The new state variable F q can be used to find the coupling with the dynamic equation of the rotor mass m . In this way it is possible to study the rotordynamic implications of supporting rotors with different types of electrodynamic bearings.

The simplest model that can be used to study the dynamic behavior of a rotor is the Jeffcott rotor model. It consists of a point mass attached to a massless shaft. This model represents an oversimplification as it neglects many aspects present in real world rotors, but, nevertheless it allows gaining insight into important phenomena especially in the case of rotors supported by EDBs.

Figure 3.

Schematic model of a Jeffcott rotor on an electrodynamic bearing showing the fundamental variables used to describe the dynamics of the whole suspension.

In this section we will study the stability of the Jeffcott rotor model supported exclusively by EDBs. The stability of a linear system is determined by its eigenvalues. Briefly, a system is stable if the real part of all the eigenvalues is negative [7]. This means that the system will exhibit a bounded output for respective bounded inputs. In the rotordynamics context this means that the rotor will respond to any disturbance forces with orbits of bounded radius.

Graphical representations are used to demonstrate the concepts, and the values of the parameters used to obtain the graphs are given in Table 1. A simplified model of a Jeffcott rotor is shown in Figure 3.

Parameter Symbol Value Unit
Rotor’s mass m r 2.025 kg
Flux linkage constant Λ 0 10 Wb/m
Bearing’s resistance R 0.286 ohm
Bearing’s inductance L 0.33 mH

Table 1.

Parameters describing the dynamics of a Jeffcott rotor on EDBs.

3.1. Undamped Jeffcott rotor

The equation of motion of the Jeffcott rotor supported by EDBs is

m r q ¨ c + F q = F ext E7

where F q is the force introduced in the system by the electrodynamic bearing and F ext is a generic disturbance force acting on the rotor’s mass. The external force can be due to gravity or rotor’s unbalance for example. Since the equation is linear, the response of the system is given by the superposition of the solution for each force. The response to a constant force, e. g. rotor weight, causes a constant displacement between the rotation axis and the symmetry axis of the magnetic field. The response to unbalance, on the other hand,is given by a whirling of the rotor and depends on the number of pole pairs [9]. Note that the force introduced by the bearing is seen as a reaction by the rotor mass.

The EDB of Eq. (6) and the rotor of Eq. (7) are interacting subsystems. The rotor responds to forces and moments with velocities and displacements. The bearing responds to the rotor’s outputs with forces. As a consequence, to study the dynamic behavior of the rotor running on EDBs, Eq. (6) and Eq. (7) must be solved together. Given the linear time invariant form of the equations a state-space model can be used for this purpose. The state space model has the form:

{ q ¨ c q ˙ c F ˙ q } = A { q ˙ c q c F q } + B { F ext } E8

The dynamic matrix A is

A = [ 0 0 1 m r 1 0 0 Λ 0 2 L j Λ 0 2 L ( p 1 ) Ω ( R L + j ( p 1 ) Ω ) ] E9

And the input gain matrix B is equal to

B= [ 1 m r 0 0 ] E10

The state-space modeling allows studying the rotordynamic stability, frequency response, unbalance response, and enables developing other tools to study the dynamics of the suspension in a fast and easy way. The analysis of different systems can be performed as simple parametric studies.

To study the rotor’s stability we calculate the eigenvalues of the dynamic matrix A of the suspension’s model (rotor supported by EDB) and analyze the evolution of the system’s poles in a root loci plot. Figure 5a shows the root loci plot obtained by calculating the eigenvalues of Eq. (9) for increasing values of rotating speed Ω . Note that the figure shows the evolution of the poles for the homopolar case and for the heteropolar with p = 2 .

It can be seen how the system presents a root that is in the right half plane for any value of rotating speed different from zero. This is true for both homopolar and heteropolar cases, representing that the Jeffcott rotor supported by EDBs is unstable for any value of rotating speed if the system is not modified. The reason for this unstable behavior has been identified to be the presence of rotating damping in the system. The eddy currents induced in the conducting part of the EDB dissipate energy associated to the motion of the rotating part. Rotating damping forces are known to destabilize the free whirling motion of the rotors for speeds above the first critical. In particular, if the rotating damping is of viscous type, the instability threshold of the undamped system (no external non rotating damping) is equal to the first critical speed [6].

Intuitively one can think that the instability arises from the fact that the system is always operating in supercritical regime because the electrodynamic supports are unable to give radial stiffness at zero rotating speed. Actually this statement is only partially valid since the behavior of the EDB cannot be correctly represented by a rotating viscous damper. The frequency dependence of the bearing’s forces must be taken into account, modifying the overall behavior. In the next sections the suspension model will be used to study the dynamic response and analyze different stabilization techniques proposed previously in the literature [1], [2], [3].

Figure 4.

Scheme of the Jeffcott rotor model with electromagnetic damping associated to rotor’s translational velocity.

3.2. Damped Jeffcott rotor

The most straightforward way to introduce non-rotating damping in the system is to do it by means of an electromagnetic damper, associating non rotating damping to the rotor’s translational velocity q ˙ c . This stabilization technique is shown in Figure 4, and has been proposed almost from the beginning of the interest in electrodynamic suspension of rotors. From the modeling point of view it consists simply in introducing a viscous damping element associated to the translational speed.

Notice that in this case the viscous damper is used as an approximated representation of the behavior of the electromagnetic damper.

As a result of this operation non-rotating damping is introduced in the model of Eq. (7) and the new equation of motion of the rotor’s mass is:

m r q ¨ c + c q ˙ c + F q = F ext E11

The dynamic matrix of the state-space model is also updated

A = [ c m r 0 1 m r 1 0 0 Λ 0 2 L j Λ 0 2 L ( p 1 ) Ω ( R L + j ( p 1 ) Ω ) ] E12

Figure 5b shows the influence of the non-rotating damping on the system’s poles. It is readily seen that the presence of damping allows stabilizing the dynamic behavior above a certain value of rotating speed Ω S . This value represents a stability threshold, being the system unstable for spin speeds below it and stable for speeds above it. Another conclusion that arises from this diagram is that the bearing’s rotating damping contribution reduces for higher values of spin speed, when the stabilizing stiffness contribution becomes dominant.

Notice that if the rotor’s spin speed is not constant it is necessary to introduce a further equation to express the dependence between angular displacement and driving torque. Since this additional degree of freedom is related to the rotation about the rotor’s axis, the rotor’s polar moment of inertia cannot be neglected. However for the present study it is acceptable to neglect this behavior and develop the study considering only constant rotational speed, thus elimination this further degree of freedom.


4. EDB’s stator on elastic supports

An alternative to the previous solution that allows introducing non-rotating damping in an effective way is to introduce a stabilizing element between the stator of the EDB and a rigid base. This element can be devised in different ways, for example, using viscoelastic elements, spring elements associated to passive eddy current dampers or even using active dampers [8]. In general, the introduction of stiffness and damping contemporarily is needed.

Figure 5.

(a) Root loci of the undamped Jeffcott rotor on EDBs; (b) root loci of the damped Jeffcott rotor on EDBs. Point markers ( · ) represent the roots for homopolar bearings; circular markers ( ) represent the roots for heteropolar bearings with two pole pairs ( p = 2 ).

Within the EDB’s context this system has been analyzed by Tonoli et al. [3]. It was shown that the stability boundaries of a Jeffcott rotor on this type of support allow operating at reduced rotational speeds with respect to the most common electromagnetic damping system proposed in literature [1], [2]. Furthermore, this stabilization technique avoids increasing the rotor’s mass and complexity because all the additional subsystems are placed on the stator part. In addition, the possibility of introducing non-rotating damping between two stationary parts allows using classical damping technologies, such as, viscoelastic materials or squeeze film dampers. However, the choice of appropriate values of stiffness and damping of the stabilizing element is not obvious, requiring the solution of an optimization problem.

Figure 6.

Schematic representation of a Jeffcott rotor supported by EDBs having elastic connections between EDB’s stator and casing of the machine.

From the modeling point of view this case is interpreted as shown inFigure 6. In the figure the x y frame is inertial as if attached to an infinitely stiff base while x ' y ' is a translating reference frame attached to the center of the stator. The axes of the two systems remain parallel to each other. It must be noticed that in this case there are three interacting subsystems, namely, the rotor, the EDB, and the stator. With respect to the figure the associated degrees of freedom are:

  • The displacement of the rotor geometric centerC in the inertial frame of reference.

q = x + j y E13
  • The displacement of the stator mass m s represented by the point S in the inertial frame.

q s = x s + j y s E14
  • The relative displacement between the displacement between point C and S.

q c = q q s E15

The equations of the rotor mass and EDB are given by Eq. (7) and Eq. (6) respectively. The displacements and speeds considered in the EDB’s equations are the relative ones ( q c and q ˙ c ). The stator’s mass dynamics is described by the following equation:

m s q ¨ s + c s q ˙ s + k s q s F q = 0 E16

The presence of the negative sign on the bearing’s force F q means that the stator mass sees this force as an external force while the rotor mass perceives it as a reaction force.

The state space model of this system can be written as:

{ q ¨ q ¨ s q ˙ q ˙ s F ˙ q } = A { q ˙ q ˙ s q q s F q } + B { F ext } E17

The dynamic matrix A of this state-space model is:

[ 0 0 0 0 1 m r 0 c s m s 0 k s m s 1 m s 1 0 0 0 0 0 1 0 0 0 Λ 0 2 L Λ 0 2 L j Λ 0 2 L ( p 1 ) Ω j Λ 0 2 L ( p 1 ) Ω ( R L + j ( p 1 ) Ω ) ] E18

And the input gain matrix B is equal to

B = [ 1 m r 0 0 0 0 ] E19

The root loci of this system considering the same bearing’s characteristics of the previous case are shown in Figure 7a. The values of stiffness k s and damping c s are 240 kN/m and 510 N s/m. The choice of these values must be done performing an optimization to minimize the stabilization threshold speed. To have an objective view of the problem one can resort to a plot showing the value of the stabilization threshold speed in terms of different values of k s and c s . Figure 7b shows a contour plot of the stabilization threshold speed. The presence of a minimum is clear in the figure leaving to the designer the task of optimizing the system’s properties to minimize the stabilization threshold speed guaranteeing that the region surrounding the minimum lies within a region of physically feasible property values.

Figure 7.

(a) Root loci of the Jeffcott rotor on EDBs having elastic elements in between the stator and the casing of the machine. (b) Contour map of the stabilization threshold speed for different sets of values of damping c s and stiffness k s of the elastic element connecting EDB’s stator and casing. The lines crossing the graph from lower left to upper right denote the value of loss factor η = c s / k s m s associated to that configuration.

4.1. Anisotropy of heteropolar bearings

In the preceding sections both rotor and stator were assumed to be axial symmetric. Considering the difficulty in insuring stability of the whirling motion of the rotor, a stabilizing technique for transverse whirl modes introducing anisotropy into the bearing stiffness can be considered [2]. This can be achieved in different ways, but one simple strategy is the use of an anisotropic Halbach array of magnets, where the gradient of the flux density in one direction is different from the other, thus modifying the value of the parameter Λ 0 in one direction relative to the other. Another possible strategy is to use rotating magnets and fixed conductors, and to have a different set of properties of the electrical circuit in each direction. For the present study only the first strategy is considered. To study this type of physical problem the state-space model given in terms of complex coordinates in Eq. (17) must be split into its representation in terms of real coordinates as:

{ x ¨ c y ¨ c x ˙ c y ˙ c F ˙ x F ˙ y } = A { x ˙ c y ˙ c x c y c F x F y } + B { F ext } E20

where the dynamic matrix assumes the form:

[ c m r 0 0 0 1 m r 0 0 c m r 0 0 0 1 m r 1 0 0 0 0 0 0 1 0 0 0 0 Λ 0 x 2 L x 0 0 ( p 1 ) Λ 0 x 2 L x Ω R x L x ( p 1 ) Ω 0 Λ 0 y 2 L y ( p 1 ) Λ 0 x 2 L x Ω 0 ( p 1 ) Ω R y L y ] E21

Calculating the eigenvalues of the dynamic matrix for different values of spin speeds it is possible to fine the stabilization threshold speed. If different values of the ration between the properties in x direction and those in y direction, and finding the stabilization threshold speed in every case it is possible to study how the anisotropy of these properties influence the stabilization speed. Figure 8 shows the graphs obtained performing this operation for different values of non-rotating damping between rotor and stator. It can be noticed that the anisotropy has a strong influence on the stabilization speed. In fact, one of the worst cases is precisely when the properties of the bearing are isotropic; this is evidenced by the peak in the stabilization threshold speed. The non-rotating damping has the effect of reducing the stabilization speed in the entire speed range. One important aspect is the evidence that for larger values of anisotropy combined to large values of bearing’s stiffness, the stabilization threshold reduces strongly towards zero also for low values of non-rotating damping.

Figure 8.

Stabilization speed of the rotor on heteropolar EDB with anisotropic properties of the bearing.

4.2. Anisotropy of stator-casing connections

The homopolar concept was first devised to eliminate unnecessary eddy-current losses generated by AC electrodynamic bearings [5]. The concept itself presupposes axial symmetry of both rotor and stator; hence the introduction of anisotropy of the bearing is not possible. On the other hand, considering the configuration presented in Sec. 4, it is possible to imagine a system where the stiffness and damping of the connection between EDB’s stator and casing are different in each direction.

Similarly to the previous case, this system is more conveniently represented in real coordinates. The representation in complex coordinates is possible as well but creates difficulties for the state-space modeling.

In the first paragraph the homopolar concept was cited to motivate this section, however, as a consequence of the unified modeling, the effect of anisotropy can be appreciated in both homopolar and heteropolar configurations. The state-space model can be written as:

{ x ¨ c y ¨ c x ¨ y ¨ x ˙ c y ˙ c x ˙ y ˙ F ˙ x F ˙ y } = A { x ˙ c y ˙ c x ˙ y ˙ x c y c x y F x F y } + B { F ext } E22

where the dynamic matrix is:

A = [ 0 0 0 0 0 0 0 0 1 m r 0 0 0 0 0 0 0 0 0 0 1 m r 0 0 c x m s 0 0 0 k x m s 0 1 m s 0 0 0 0 c y m s 0 0 0 k y m s 0 1 m s 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 Λ 0 2 L 0 Λ 0 2 L 0 0 Λ 0 2 L Ω 0 Λ 0 2 L Ω R L Ω 0 Λ 0 2 L 0 Λ 0 2 L Λ 0 2 L Ω 0 Λ 0 2 L Ω 0 Ω R L ] E23

From the stability point of view the inputs of the linear system are irrelevant and the input gain matrix doesn’t have to be defined.

Figure 9.

Stabilization speed of the rotor on homopolar EDB with anisotropic connections between bearing stator and casing of the machine.

To study the possibility of taking advantage of anisotropy of the connections to reduce the stabilization threshold speed, the stabilization speed is calculated for different values of the anisotropy ratio ( α = k x / k y ). Considering constant values of rotor and stator masses, m and m s respectively, and the system scheme of Figure 6, the effect of anisotropy is studied in terms of the damping ratio ζ . Figure 9 shows how anisotropy and damping ratios affect the stabilization speed. For systems having low damping ratios the anisotropy can have a beneficial role, reducing the stabilization threshold speed.Increasing the damping ratio eliminates the positive effects due to anisotropy in the elastic connections, but effectively reduces the stabilization threshold.

The anisotropy in this case has a different effect with respect to that illustrated in Figure 8. The increase in the anisotropy ratio with a respective increase of stiffness of one direction results to increase the stabilization threshold speed. It is obvious that this diagram is case dependent, but in general it is expected that the anisotropy has a positive contribution only when the value of damping is low [6]. Furthermore, within physically feasible margins it is always more advantageous to increase the value of damping than to use effects of anisotropy because the stabilization threshold is more sensitive to the first than to the latter.


5. Conclusions

The present paper presents the development of a dynamic model of the radial suspension using electrodynamic bearings that is adopted to study the mechanical properties of the supports that allow guaranteeing rotordynamic stability. A simple procedure is used to identify the characteristics of the bearing, in case of heteropolar bearings, and of the elastic support that allow obtaining the best performance in terms of minimization of stabilization speed.

The effect of anisotropy of the supports in the stabilization threshold speed is also investigated. It is noticed that the anisotropy of the EDB’s properties in case of heteropolar configurations can be advantageous independently of the amount of nonrotating damping that can be introduced. The anisotropy allows obtaining stabilization speeds that are lower than the isotropic case. In fact the isotropic bearing represents a critical case, with extremely high stabilization speeds with respect to an anisotropic configuration.

In case of homopolar EDB configurations it is not possible to devise an anisotropic bearing because of the intrinsically axisymmetric distribution of the magnetic field. Hence the anisotropy of the elastic elements connecting the EDB’s stator to the casing of the machine has been analyzed under the same hypothesis assumed in the case of heteropolar configurations. It has been observed that anisotropic characteristics of the supports can be advantageous only at low damping levels. For higher values of damping of the connection element the advantages of anisotropy vanish, and the isotropic configuration becomes optimal. Furthermore, it has been observed that is more advantageous to increase damping instead of resorting to anisotropic configurations in the case of anisotropy of the elastic connections because the stabilization threshold speed is more sensitive to the first than to the latter.


  1. 1. Post R. F. Ryutov D. D. 1998 Ambient-temperature passive magnetic bearings: Theory and design equations Proceedings of the 6th International Symposium on Magnetic Bearings Cambrige Massachusetts, USA
  2. 2. Filatov A. Maslen E. H. 2001 Passive magnetic bearing for flywheel energy storage systems IEEE Transactions on Magnetics 37 6
  3. 3. Tonoli A. Amati N. Impinna F. Detoni J. G. 2011 A solution for the stabilization of electrodynamic bearings: modeling and experimental validation ASME Journal of Vibration and Acoustics 133
  4. 4. Lembke T. A. 2005 Design and analysis of a novel low loss homopolar electrodynamic bearing PhD thesis Royal Institute of Technology Stockholm, Sweden
  5. 5. Kluyskens V. Dehez B. 2009 Parameterized electromechanical model for magnetic bearings with induced currents Journal of System Design and Dynamics 3 4
  6. 6. Genta G. 2005 Dynamics of rotating systems Springer New York
  7. 7. Dorf R. C. Bishop R. H. 2010 Modern Control Systems Prentice Hall
  8. 8. Tonoli A. Amati N. Bonfitto A. Silvagni M. Staples B. Karpenko E. 2010 Design of Electromagnetic Dampers for Aero-Engine Applications Journal of Engineering for Gas Turbines and Power 132 11
  9. 9. Detoni J. G. Impinna F. Tonoli A. Amati N. 2012 Unified Modelling of Passive Homopolar and Heteropolar Electrodynamic Bearings Journal of Sound and Vibration 331 19 4219 4232

Written By

J. G. Detoni, F. Impinna, N. Amati and A. Tonoli

Submitted: March 8th, 2012 Published: October 2nd, 2012