Self-Identification Algorithm for the Autonomous Control of Lateral Vibration in Flexible Rotors

Intelligent machines are capable of recognizing their operational condition and take actions towards optimality through an autonomous processing of information. Considering the importance of rotating machines in modern industry, this concept of intelligent machines can be applied to achieve high availability, thus avoiding interruptions in the production flow. In this work, a self-identification algorithm is proposed for the autonomous decision and control of a flexible shaft rotating system with electromagnetic actuators. Based on the D-decomposition technique, the algorithm searches in the domain of controller gains the best ones for P and PD controllers to reduce maximum peak response of the shaft. For that, frequency response functions of the system are automatically identified experimentally by the algorithm. It is demonstrated that regions of stable gains can be easily plotted, and the most suitable gains can be found to minimize the resonant peak of the system in an autonomous way, without human intervention.


Introduction
According to [1], an intelligent machine is capable of recognizing its operational condition and takes actions towards optimality through an internal and autonomous processing of information.This can be understood as a system able to adapt and learn from plant data and assure the desired dynamic behavior of the plant under any disturbance and/or change in its inner characteristics [2].
In the last decades, actuators and sensors have been incorporated into rotating machines aiming at attenuating and controlling the vibration levels, especially those presented by the shaft [3][4][5][6].In such cases, the actuators can also be used as exciters, and open-loop frequency response of the system can be identified.Such frequency response would contain information from the actuator system + rotating system/ bearings + sensor system that can be used to find the best gains for the controller.Therefore, if all this process (identification of plant FRFs and determination of optimum gains) is automatically managed by an algorithm, as proposed in the present work, one would achieve the definition of intelligent machines stated by [1]: a machine that can identify its characteristics (self-identification) and take actions towards optimality (search for optimum gains).
In the area of active magnetic bearings (AMBs), much has been done in this direction by the development of adaptive control systems.In most cases, the conventional control system, responsible for rotor levitation, is complemented by an additional adaptive algorithm, responsible for adjusting control parameters to external disturbances.Examples of such procedure are found in [7][8][9], where the adaptive controller is designed to attenuate rotor unbalance, which is usually an unknown disturbance.More recently, the work of [10] presents a multiobjective adaptive controller, responsible for reducing vibration in different parts of the machine as well as reducing the transmitted force to the foundation.In all these works, the conventional control system is designed previously, and the rotating system adapts itself to new disturbing conditions during operation.
In this work, the idea is not adapting to new environment conditions, but finding the best gains of the conventional controller in an autonomous way, with no previous International Journal of Rotating Machinery information.The control system (actuators and sensors) shall identify the system characteristics (dynamics) and decide which gains to use in a PD controller (conventional), autonomously.After this, the gains remain constant, but the procedure could be reused to update the gains after some time interval.Alternatively, the procedure could be used to find the first gains to be adopted in the conventional controller and, later, the gains could be updated by an adaptive controller according to the external disturbances.Considering that AMBs need a controller to operate from the beginning (to levitate and start running the rotor), the proposed procedure is more suitable for systems that can start running the rotor without necessarily turning on the controller: active squeeze film dampers [11], bearings with active lubrication [5], bearings with piezo stack actuators [6], active hydrostatic bearings [12], and any other kind of active hybrid bearing.
One of the first works on controller synthesis based on frequency response functions of the plant refers to linear controllers [13].In this case, the linear controller is similar to a lead/lag compensator, and the closed-loop stabilizing gains of the controller are determined from the frequency response functions of the open-loop system.Later, the theory for designing PID controllers was presented [14], whose equations for finding the stabilizing gains were formally deduced in [15].By knowing the gains that stabilize the closed-loop system, one can find sets of gains that optimize any performance criteria, thus achieving closed-loop stability and robustness [16].
The cases cited above adopt the D-decomposition technique to find the regions of stabilizing gains.The technique consists in deriving three conditions from which one can decompose the parameter space into regions with a fixed number of stable and unstable roots (root invariant regions).Recently, the D-decomposition technique has been used together with mathematical models to find the gains of linear and H ∞ controllers [17][18][19].However, few works in the literature apply the D-decomposition technique to find the stabilizing gains from measured frequency response functions (FRFs).
The application of such methodology to rotating systems is a novelty and it is still open for investigations.A pioneer work in this sense is presented in [20], where such model-free control technique is experimentally applied to a Laval rotor, in order to ensure stability and performance of the closedloop system.In this work, the technique is implemented in an algorithm for the autonomous identification and controller synthesis of a rotating machine.In this algorithm, the Ddecomposition technique is used to find the stabilizing gains of P and PD controllers, whose objective is controlling lateral vibration of the shaft.The methodology is implemented experimentally and regions of stabilizing gains are obtained from the global frequency response of the openloop system (actuators + rotor/bearings + sensors).The optimum gains are selected to attenuate vibrations according to a performance criterion, and implemented automatically.The obtained results show the feasibility of the algorithm in finding the best gains for the control system in an autonomous way.

Performance criterion
Search for optimum gains

Implementation of optimum gains
Calculation of optimum gains Figure 1: Flow chart of the algorithm for autonomous selfidentification and control of the rotating system.

The Self-Identification Algorithm
The idea of the self-identification algorithm is to employ the actuation system to identify the system dynamics and, from this information, find the appropriated gains to control the system.A general flowchart of the algorithm is shown in Figure 1.Initially, the controller is turned off (null gains) and a subroutine starts the identification procedure, which consists in measuring the frequency response function (FRF) of the open-loop system (actuator + rotor/bearings + sensors).Next, the D-decomposition technique is used to find the gains that stabilize the system based on the measured FRF of the system.By adopting a performance criterion for controlling shaft vibration (e.g., minimization of resonance peak amplitude, minimization of control power, minimization of averaged vibration levels along the shaft), the algorithm searches for optimum gains among those defined by the Ddecomposition technique.Finally, the controller is turned on with the optimum gains found.All procedure is done autonomously by the machine (algorithm).

The D-Decomposition Technique Applied to PD Controllers
The basic idea of the D-decomposition technique consists in dividing the controller parameter space into root invariant regions, that is, regions with a fixed number of stable and unstable roots [21].This division is obtained calculating boundaries which map the imaginary axes in the complex plane into curves in the controller parameter space.Because each region delimited by these boundaries corresponds to a root invariant region, if a set of gains in a certain region leads the system to stability, all gains inside this region will do it as well.
Consider π( jω, λ) as the characteristic polynomial of a closed-loop system of degree n written in the frequency domain and with real coefficients α k (λ), where λ ∈ m is an uncertain parameter (e.g., controller gains) and ω ∈ [0, +∞).The space formed by the controller gains can be divided into regions denoted by D(m, n − m), which correspond to the characteristic polynomial with m roots with negative real part (stable poles) and n − m roots with positive real part (unstable poles).This decomposition in the parameter space into regions D(m, n − m) is called D-decomposition.The boundaries of each region are defined by the poles in the stability margins of the system (Figure 2), given by Equation ( 1) is known as the real root boundary (RRB) and provides the values of the parameter λ for which the real roots are in the origin.Equivalently, it can be written as Equation ( 2) is known as the complex root boundary (CRB) and provides the values of the parameter λ for which the complex roots are on the imaginary axis.Equation ( 3) is known as the infinite root boundary (IRB) and provides the values of the parameter λ for which the roots pass through infinity, that is, the other way the real roots can cross the imaginary axis (through a reduction in the order of the characteristic polynomial thus consequently changing the number of the poles).
The three previous conditions allow one to decompose the controller parameter space in regions with a fixed number of stable and unstable roots (root invariant regions) [21].However, there is no guarantee of the existence of a region D(n, 0), that is, a region with stable roots only, and all the obtained regions can result in closed-loop system instability [22].

PD Controller Root Invariant Regions.
Consider the feedback configuration with a linear time invariant system and a first-order controller according to Figure 3.The system and the controller are defined, respectively, as where G P and G D are the proportional and derivative gains of the PD controller, and P( jω) is the frequency response of the system, which is assumed to be controllable.The closed-loop response function of the system in the frequency domain is given by Thus, the characteristic polynomial of the closed-loop system is where subscripts r and i refer to the real and imaginary parts, respectively, and Hence, from ( 7) and ( 8), ( 1) to (4) can be written as functions of the frequency response of the open-loop system and of the gains of the controller.
(a) The Real Root Boundary (RRB).Considering the closed loop characteristic polynomial of the system (7) and the real root boundary defined by ( 1) and ( 4), one has If P(0) / = 0, then

(b) The Infinity Root Boundary (IRB).
Assuming that the denominator of the plant's transfer function has order n and the numerator has order n − 1, and taking the infinity root boundary defined by (3), one can write the closed loop characteristic polynomial (7) as where d n and n n are, respectively, the nth-order coefficients of the denominator and the numerator.Hence, if P(∞) / = 0:

Test Apparatus
In the case P( jω) is obtained experimentally, the infinite root boundary is seldom applicable because either information at infinity is not available or system response tends to zero (G D −→ ∞).
(c) The Complex Root Boundary (CRB).According to the complex root boundary defined by ( 2), ( 8) can be written as whose solution is Hence, given the frequency response of the system to be controlled P( jω), the PD controller parameter space (G D × G P domain) can be partitioned into root invariant regions whose boundaries are defined by ( 10), (12), and ( 14).The only information required for the analysis is the frequency response of the open-loop system.In order to know which partitions correspond to D(n, 0) regions, that is, regions with n stable and zero unstable poles, it is necessary to select a point from each region and test the stability of the corresponding closed-loop system.
The proposed method is experimentally applied to the double disk rotating system shown in Figure 4.The test rig is composed of a steel shaft, supported by self-aligning ball bearings, with two steel disks located at ∼1/6 of the shaft span between the bearings (gyroscopic effects are expected).The shaft is driven by an electric motor whose speed is controlled by a frequency inverter.Electromagnetic actuators are located near the disk no.1, acting in the horizontal direction (perpendicular to the view shown in Figure 4).The response of the shaft is measured by proximity probes mounted near the two disks in the horizontal and vertical directions, but only measurements in horizontal direction near disk #1 are used for control system feedback.A NI PCI-6229 acquisition board is used for signal acquisition and control.
To obtain the frequency response of the system, the actuating system + shaft dynamics + sensor system will be considered as a whole, whose response is P( jω) (Figure 5).In this case, all dynamics of the actuating system (drives and actuators), rotating system and sensor system will be embedded in the measured frequency response P( jω).The open-loop frequency response of the system is presented in Figure 6.The fundamental frequency of the shaft associated to the first bending mode is ∼26 Hz, and the resonance associated to the second bending mode is beyond the frequency range of study (above 70 Hz).Gyroscopic effect is responsible for the separation of the resonance in two peaks as rotating speed increases.The 1 × Ω component in Figure 6 is caused by residual unbalance.

Example of an Implemented Algorithm
The self-identification algorithm is implemented considering that it is possible to find the stabilizing gains of PD controllers from frequency response information (no mathematical model required).One possibility of implementing the algorithm is shown in Figure 7. Initially, the controller gains are set to zero and the counter variable aux is initiated.The first step is identifying the system dynamic response.This first information about the system is obtained by subroutine FRF Measurement, where the frequency response functions (FRFs) of the system are identified for null controller gains (open-loop response).
The next step is finding the best stabilizing gains to be used in the closed-loop control system.For that, in subroutine Gain Search, the limits for stable controller gains (G P , G D ) are found according to the D-decomposition technique ((10) to ( 14)).The performance criterion of resonance peak amplitude is chosen, and a search algorithm is implemented to find the gains that minimize this criterion.
After finding the gains that minimize the resonance peak amplitude of the closed-loop system, the gains are implemented.For that, every pair of gains (G P , G D ) is successively implemented starting from (G P = 0, G D = 0), and the FRF of the closed-loop system is identified.The identification of closed-loop FRFs is important to check if the control signal does not saturate the D/A port of the acquisition system.The adopted acquisition board (NI PCI 6229) limits the output signals at the A/D and D/A ports to ±10 V.
When the D/A port is saturated by the control signal (u c ), it means that the control gain values are too high for the application.Hence, a strategy of fine-tuning of the gains is adopted: the next gains to be tested result from the mean value of the last two gains already implemented.If this new gains result in a control signal that does not saturate the D/A port, then these are the final gains to be adopted in the closed-loop system.If not, the final gains to be adopted in the closed-loop system are those that did not saturate the D/A port one iteration before (aux-1).

FRF Measurement Subroutine. The sub-routine FRF
Measurement is responsible for the identification of system FRFs (Figure 8).For that, a chirp signal is sent to the actuators and system response is measured by the proximity probes.The chirp signal varies from 1 to 70 Hz in 5 s and returns to 1 Hz in another 5 s.This signal is repeated four times, resulting in a 40 s period of data acquisition.The adopted acquisition rate is 1 kHz.
After measuring the input (chirp signal) and system output (displacements from the proximity probes), the FRFs are calculated using the H 1 estimator [23].FRF (Figure 9).Initially, the invariant root regions in the G P × G D domain are calculated, from (10) to (14).Equations ( 10) to ( 14) result in curves in the G P × G D domain that delimit the range of gains that stabilize the system (limiting values for G P and G D before the closed-loop system becomes unstable).

Gain Search
Once the region of stable roots is determined, one can choose a set of gains inside this region which satisfies a desired criterion of performance.For instance, consider the resonance peak amplitude as such criterion.Based on the idea of using only the frequency response of the system in the controller design (model free approach), consider the closedloop transfer function of the system in the frequency domain given by ( 6) Varying the values of G P and G D , one can calculate the expected resonance peak amplitude of the closed-loop system for the gains of the controller: where the resonant peak to be minimized is within the range of frequencies defined by ω max and ω min .The resonance peak amplitude of the closed-loop system (M r ) is then calculated (16) within the limiting range of G P and G D for a stable closed-loop system.This gives a surface in G P × G D domain that can be used to find the gains that minimize M r .The gains that minimize M r are found in windows of Figure 7: Flow chart of the implemented algorithm for autonomous self-identification and control of the rotating system.
forming a sequence of gains towards the minimum of M r surface, starting from G P = 0, G D = 0 (Figure 10).After 11 iterations (windows), the calculation of the best gains stops.The maximum number of iterations depends on the adopted value for the step variable and the limits of G P and G D for stable closed-loop system.Mathematically, the sequence of gains represents the maximum negative gradient of surface M r .Hence, one can infer that this procedure could be done without the windows.However, surface M r comes from experimental data and, as a consequence of the inherent noise in the measurements, it is not a smooth surface.Experience demonstrates that finding the sequence of gains from the calculated gradient of surface M r does not achieve good results, and this procedure shall only be used in the case of low noise measurements.

GUI Interface.
The algorithm was implemented and a user interface was developed with the GUI tools of MatLab (Figure 11).By clicking on the start button, the program runs

Experimental Results
The self-identification algorithm was implemented and tested for two different rotating speeds of the test rig: subcritical speed of 280 rpm (4.7 Hz) and supercritical speed of 2080 rpm (34.7 Hz).In Figure 13, one presents the obtained frequency response functions for the subcritical rotating speed of 280 rpm at sensor #1 (lateral vibration of the shaft).The maximum peak amplitude of the open-loop system (black line) has a value of −24.85 dB (0.057 mm/V) (dB = 20 log 10 (mm/V)) at ∼26 Hz.The gyroscopic effect at this rotating speed is weak, and the resultant separation of the fundamental frequency is not evident (second peak is small).By implementing the PD controller (dark gray line), the algorithm found gains that reduced the resonance peak amplitude in 27% to the value of −27.61 dB (0.042 mm/V) and shifted the resonance peak to ∼27 Hz.Another noticeable effect of the PD controller is that it increased the amplitude of the second peak of the resonance, and one can clearly see the separation of the fundamental frequency due to the gyroscopic effect (two peaks).A third attempt of the algorithm was based on the search of gains for a P controller (only proportional gains were searched and implementedlight gray line).As a result, the resonance peak amplitude reduced 54% to the value of −31.57dB (0.026 mm/V) and it is not observed a shift in resonance frequency or any amplification of the second peak.
The gains autonomously tested by the algorithm are presented in Figure 14.In this figure, the numbers represent the sequence of tested gains (obtained by sub-routine Gain Search).For the PD controller, the adopted step for the windows was 0.025, which represents incremental steps of 0.025 Vs/mm for the derivative gain G D and 0.5 V/mm for the proportional gain G P .For the P controller, the adopted step for the windows was 0.1, which represents incremental steps of 2.0 V/mm for the proportional gain G P .
The choice of the incremental step in the gain search procedure is a trade-off between the desired maximum number of iterations and the appropriated tuning of the final gains.If the incremental step is too small, an excessive number of gains will be tested by the algorithm until it finds the final gains, but the gains will be finely tuned.If the incremental step is too big, the number of iterations is strongly reduced, but fine-tuning of the gains will be jeopardized (difference between tested gain values will be large).In the present work, the chosen incremental steps resulted in a small number of iterations with a reasonable fine-tuning of the gains (Figure 14).
For the PD controller (dark gray line), the algorithm tested the best gains until the fifth iteration, when there was saturation of the D/A port.For this reason, the algorithm tested values for the gains between those tested in the fourth and the fifth iterations (point 6), which still saturated the D/A port.As a result, the final gains adopted by the algorithm were those tested in the fourth iteration (point 7 ≡ point 4).For the P controller (light gray line), the algorithm tested the gains until the fourth iteration, when there was saturation of the D/A port.The final adopted gains were those obtained by the fine-tuning: the mean value between the third and fourth gains tested.
By comparing the results shown in Figure 13, the better performance of the P controller is explained by the fact that the algorithm "walked farther" in the G P × G D domain towards the minimization of M r in comparison to the sequence of the PD gains (Figure 14).This happened because the PD controller needs shaft velocity information, and derivation of displacement signals from the proximity probes introduced high noise level in the control signal.As a result, the D/A port was more easily saturated by the PD control signal than by the P control signal, and higher values of gains could be used in the P controller case.An improvement in the signal-to-noise ratio of shaft velocity data will surely allow the algorithm to test gains beyond point 5 in Figure 14, thus reducing the amplitude of the closed-loop resonance peak with the PD controller.
The invariant root regions of the system for the rotating speed of 280 rpm are shown in Figure 15.As one can see, two curves delimit three areas in the G P × G D domain.Considering that the rotating system supported by ball bearings is inherently stable, the inner region which contains the origin is a region of stable roots (this can also be proven by a stability analysis).The shaded area in Figure 15 represents the area of tested gains by the algorithm, which remained within the limits of the stable region.Further calculation of optimum gains is possible until the limits of the windows reach the borders of the invariant root regions in the G P × G D domain.However, saturation of D/A ports in the acquisition system limits the range of tested gains.The only data used to feed back the control system are those of sensor #1 in horizontal direction.Despite this unidirectional control, the frequency response functions of the system in the vertical direction were not strongly neither negatively affected (Figure 16).
In Figure 17, one presents the obtained frequency response functions for the supercritical rotating speed of 2080 rpm at sensor #1.At this rotating speed, the gyroscopic effect is stronger, and the resultant separation of the fundamental frequency is clearly noticeable (two peaks related to the same mode shape: first bending mode).The maximum amplitude of the peaks for the open-loop system (black line) has values of −29.55 dB (0.033 mm/V) at ∼23 Hz (first peak), and −28.95 dB (0.036 mm/V) at ∼28 Hz (second peak).By implementing the PD controller (dark gray line), the algorithm found gains that increased the amplitude of the first resonance peak in 9% to the value of −28.93 dB (0.036 mm/V) at ∼22 Hz and reduced the second resonance peak in 22% to the value of −31.15 dB (0.028 mm/V) at ∼28 Hz.The implementation of the P controller resulted in a  reduction of the first resonance peak amplitude in 42% to the value of −34.37 dB (0.019 mm/V) at ∼23 Hz, and a reduction of the second resonance peak amplitude in 42% to the value of −33.44 dB (0.021 mm/V) at ∼28 Hz.
The component at 34.7 Hz is caused by residual unbalance.This strong 1 × Ω component appears in the FRF because the chirp excitation signal, which is a time-varying frequency signal, is not correlated to the constant frequency component in the response signal caused by the residual unbalance.However, it is expected a reduction of unbalance vibration as soon as the P or PD controllers take action.Hence, the presence of a strong 1 × Ω component in the  closed-loop FRFs represented a drawback in the identification procedure at these specific frequencies (rotating frequencies).The gains autonomously tested by the algorithm at the rotating speed of 2080 rpm are presented in Figure 18.For the PD controller (dark gray line), the algorithm tested gains until the third iteration, when there was saturation of the D/A port.The final gains adopted by the algorithm were those obtained by the fine-tuning, in between the second and third iterations.For the P controller (light gray line), the algorithm tested gains until the fifth iteration, when there was saturation of the D/A port.The final gains adopted by the algorithm were those obtained by the fine-tuning, in between the fourth and fifth iterations.
Again, the better performance of the P controller in comparison to the PD controller is explained by the fact that the algorithm "walked farther" in the G P × G D domain towards the minimization of M r in comparison to the PD controller (Figure 18).The derivation of displacement signals at higher rotating speeds introduces higher noise levels in the PD control signal than those presented at the rotating speed of 280 rpm.As a result, the sequence of tested PD gains was even shorter than that tested at the rotating speed of 280 rpm, and the effects of the PD controller on the attenuation of lateral shaft vibration were not significant.
The invariant root regions of the closed-loop system for the rotating speed of 2080 rpm is very similar to those of the system for the rotating speed of 280 rpm, shown in Figure 15, and the area of tested gains by the algorithm remained within the limits of the stable region.The only data used to feed back the control system are those of sensor #1 in horizontal direction.Despite this unidirectional control, the frequency response functions of the system in the vertical direction were not strongly neither negatively affected (Figure 19).

Conclusion
The application of the proposed self-identification algorithm, based on the D-decomposition technique, for the autonomous control of a flexible shaft rotating system with electromagnetic actuators resulted in the following conclusions.(i) The algorithm successfully identifies the system dynamic characteristics, finds the best controller gains, and implements the controller in an autonomous way.
(ii) The identification process and definition of the best gains is performed without a mathematical model of the system.In fact, all system information is contained in the frequency response function of the actuator + plant + sensor, as a whole.Therefore, modeling uncertainties and misleading assumptions are avoided in the controller gain determination, and only actual information is considered.
(iii) The P controller presented better results when compared to the results obtained with the PD controller.
International Journal of Rotating Machinery The poor performance of the PD controller can be explained by the fact that PD controllers need velocity information, and the derivation of displacement signals from the proximity probes introduces high noise level in the control signal.Hence, the control signal easily reaches the saturation level of the D/A ports of the acquisition board, and one cannot adopt higher values for the control gains.If the actuators had higher actuation capacity or the derivation of displacement signals presented a better signal-to-noise ratio, higher control forces would be achieved for a given set of controller gains, and the PD controller could present a better performance.
The methodology adopted in this work involved excitation and actuation in one single direction and the results in the perpendicular direction of actuation do not show degradation of vibration levels (probably due to low gyroscopic effects).However, the adoption of the proposed methodology in bidirectional actuation must be further investigated.

Figure 2 :
Figure 2: Poles at the stability margins of the closed-loop system.
−step G D step and −20 step G D 20 step.For every window, the gains G P , G D that minimize M r are found, thus International Journal of Rotating Machinery |u c | > lim |u c | > lim

Figure 8 :Figure 9 :
Figure 8: Flow chart of the subroutine FRF Measurement.

1 PathFigure 10 :
Figure 10: Contour plot of M r in the G P × G D domain: sequence of gains that minimize M r within windows.

Figure 11 :Figure 12 :
Figure 11: Graphical user interface of the self-identification algorithm based on MatLab.

Figure 13 :
Figure 13: Frequency response functions identified by the algorithm for the rotating speed of 280 rpm.

Figure 14 :
Figure 14: Controller gains tested by the algorithm for the rotating speed of 280 rpm and contour plot of M r .

Figure 15 :
Figure 15: Invariant root regions for the rotating speed of 280 rpm and area of tested gains.

Figure 16 :
Figure 16: Frequency response functions of the system in the vertical direction for the rotating speed of 280 rpm (sensor #1).

Figure 17 :
Figure 17: Frequency response functions identified by the algorithm for the rotating speed of 2080 rpm.

Figure 18 :
Figure 18: Controller gains tested by the algorithm for the rotating speed of 2080 rpm and contour plot of M r .

Figure 19 :
Figure 19: Frequency response functions of the system in the vertical direction for the rotating speed of 2080 rpm (sensor #1).