Next Article in Journal
A Result Regarding Finite-Time Stability for Hilfer Fractional Stochastic Differential Equations with Delay
Previous Article in Journal
Effects of Transverse Friction Massage on the Electromechanical Delay Components and Fractal Dimension of Surface Electromyography in Quadriceps Muscles
Previous Article in Special Issue
Fractional Modeling and Control of Lightweight 1 DOF Flexible Robots Robust to Sensor Disturbances and Payload Changes
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Fractional-Order Models of Damping Phenomena in a Flexible Sensing Antenna Used for Haptic Robot Navigation

by
María Isabel Haro-Olmo
1,*,
Inés Tejado
2,*,
Blas M. Vinagre
2 and
Vicente Feliu-Batlle
3
1
Instituto de Investigaciones Energéticas y Aplicaciones Industriales, Universidad de Castilla-La Mancha, 13005 Ciudad Real, Spain
2
Escuela de Ingenierías Industriales, Universidad de Extremadura, 06006 Badajoz, Spain
3
Escuela Técnica Superior de Ingeniería Industrial, Universidad de Castilla-La Mancha, 13071 Ciudad Real, Spain
*
Authors to whom correspondence should be addressed.
Submission received: 29 June 2023 / Revised: 28 July 2023 / Accepted: 29 July 2023 / Published: 15 August 2023
(This article belongs to the Special Issue Fractional Order Control Techniques for Robots)

Abstract

:
In this paper, two types of fractional-order damping are proposed for a single flexible link: internal and external friction, related to the material of the link and the environment, respectively. Considering these dampings, the Laplace transform is used to obtain the exact model of a slewing flexible link by means of the Euler–Bernoulli beam theory. The model obtained is used in a sensing antenna with the aim of accurately describing its dynamic behavior, thanks to the incorporation of the mentioned damping models. Therefore, experimental data are used to identify the damping phenomena of this system in the frequency domain. Welch’s method is employed to estimate the experimental frequency responses. To determine the best damping model for the sensing antenna, a cost function with two weighting forms is minimized for different model structures (i.e., with internal and/or external dampings of integer- and/or fractional-order), and their robustness and fitting performance are analyzed.

1. Introduction

Haptics is the technology of touch. It has found many applications in robotics such as accurate positioning, mobile robot navigation, robotic medical exploration and surgery, and collaborative robotics. These tactile sensors are complex and exhibit a distributed nature. Then, issues such as placement, sensor robustness, or wiring complexity have to be addressed. Often, sensors used in haptics try to mimic nature. In particular, haptic sensors inspired by insect antennae or mammal whiskers (see Figure 1) have been built since the early 1990s.
The “sensing antenna” has become the most popular among these sensors because it is a robust and compact device. It is an active sensor composed of a very slender flexible beam, which is moved by two servo-controlled motors, and a load cell placed between the beam and the motors. This device replicates the touch sensors that some animals possess and performs an active sensing strategy in which the servomotor system moves the beam back and forth until it hits an object. At this instant, information on the motor angles, combined with force and torque measurements, allows us to calculate the positions of the contact points, providing valuable information about the object surface. Thus, a 3D map of the surface of an object can be obtained using this device, enabling its recognition. Two strategies can be applied to obtain such 3D maps. The first one is to keep moving the beam back and forth in order to hit the object at different points, allowing for the determination of their 3D coordinates and, subsequently, extracting the map of the object surface. This procedure is carried out by some insects that use their antennae for this purpose (e.g., [1]). The other strategy is to slide the beam across the object, exerting a controlled force on the surface of the object in order to maintain contact, and collect the 3D coordinates of points on the object surface during this movement. This approach is utilized by some mammals that have whiskers as sensors (e.g., [2]). Both of these strategies can be implemented using the aforementioned sensing antennae.
We are developing a sensing antenna system to aid the navigation of the prototype of a mobile robot shown in Figure 2. Our system will localize obstacles by carrying out the first of the aforementioned strategies. We obtain the 3D information about the contact point using the method reported in [3]. This method is based on the fact that the frequencies of the vibrations of the antenna change from the free motion case to the constrained motion case and, in this last case, also change as a function of the location of the contact point. Then, the torque measurements of the vibrations produced in the beam of the antenna during its movement in a recognition task are processed to estimate the first two vibration frequencies and, thereafter, estimate the contact point.
The dynamics of these antennae is quite complicated, which hinders fast and accurate movements: beam elasticity produces residual vibrations during the antenna movement and when the target has been reached. Therefore, the design of the controllers of the motors must consider these dynamics to achieve accurate and fast approaches to the target points to be investigated. Moreover, an inefficient cancellation of the beam vibrations could produce permanent collisions with the object, in which the antenna hits back and forth. These introduce delays in the recognition, diminish the quality of the estimates of the 3D points on the object surface and, consequently, reduce the efficacy of the device.
Active control can be used to remove these beam vibrations. It allows us to search for special points in a precise manner in such a way that the maximum amount of task-relevant information is provided. In order to design such high-performance control systems that drive the beam and position its tip with high precision while preventing vibrations, very accurate models of the dynamics of the antenna in free movement are required.
This article is devoted to obtaining analytic models of the dynamics of an antenna rotating in a horizontal plane, i.e., having azimuthal movements. This system is linear but has an infinite number of nearly undamped vibration modes. In particular, modeling damping phenomena in this antenna is of the utmost importance in order to achieve faster control systems, since damping increases the relative stability of the system but slows down its response. We will show that using fractional derivatives to model the damping of the antenna enhances the accuracy of the dynamic model and, based on this, we will develop new fractional-order models that describe better the antenna’s behavior.
The remainder of the paper is organized as follows. Section 2 presents the state of the art in modeling flexible link antennae. Section 3 describes the experimental setup considered in this work, consisting of a flexible sensing antenna. Section 4 deals with the dynamic model that describes the behavior of the flexible antenna. In Section 5, the methods used to identify the dynamics of the system is detailed. Section 6 shows and discusses the results obtained after identification, highlighting the models that describe more accurately the system behavior. Finally, concluding remarks are drawn in Section 7.

2. State of the Art in Modeling Flexible Link Antennae

The dynamics of a sensing antenna is mostly characterized by the dynamics of its flexible link. Two approaches to modeling flexible links have been developed: assuming a massless link with an equivalent mass lumped at its tip, and assuming mass distributed along the entire beam. In the first approach, the dynamics includes a single vibration mode, exhibiting a minimum phase behavior. For the latter approach, the dynamics of non-collocated input-output pairs have non-minimum phase behaviors whose models are described by high-order coupled non-linear differential equations. In addition, the system is driven by small servomotors in which the limited torques, intermittent operation, and the strong non-linearity owing to the static friction constrain the motor control design. In these systems, spillover effects may easily appear, i.e., instability of the closed-loop control system caused by non-modeled high-frequency dynamics and parameter uncertainties.
Taking into account the distributed mass approach, the dynamic models of flexible links can be obtained using a wide variety of methods and considering different types of damping terms. As explained in [4], the response of a flexible link can be obtained with exact methods, i.e., modal analysis or Laplace transformations, and approximate methods, i.e., assumed modes or Galerkin’s method. Regarding damping phenomena, many varied models have been considered in other works, both linear and non-linear. When integer-order is considered, examples can be found in [5,6]. In [6], two of the damping models proposed are viscous air damping, assumed to be proportional to velocity, and Kelvin–Voigt damping, which is proportional to the fourth spatial derivative of the deflection of the link. If the flexible link behavior is extended to the fractional order, different kinds of fractional vibrations can be contemplated. In [7], six classes of fractional vibrators can be found based on what terms (inertia, friction, and/or position) are defined as fractional. Only three of these classes considered fractional friction terms explicitly. For example, the authors in [8,9] define models in which the inertia and friction terms are fractional.
Within fractional vibration models that only include fractional friction terms, two kinds of common linear damping can be considered. The first one defines Kelvin–Voigt damping as fractional-order damping [10,11,12]. In [10], the response of a Euler–Bernoulli beam under quasi-static and dynamic loads is studied. In [11], the Euler–Bernoulli theory and modal analysis are used to model a slewing flexible beam with a mass at its tip and an external load. In addition, the model of a non-linear simply supported beam with harmonic excitation is derived by the Galerkin approximation method in [12]. On the other hand, a fractional-order viscous damping (dependent only on the time derivative) is considered in [13,14,15,16,17,18]. In [13], a model of a continuous beam is obtained with modal analysis, setting the fractional derivative order of the damping to 0.5 . In [14], the Adomian decomposition method is utilized to solve a simply supported beam with a load. In [15,16], an Euler–Bernoulli beam with moving load is solved by combining the modal analysis method with the Laplace transform. In [17] and [18], the homotopy perturbation method and the homotopy analysis method, respectively, are used to obtain the dynamic response of a beam with external loads. To the authors’ knowledge, no model of a flexible link with both fractional-order dampings has been published yet. This model is studied in this work, and an exact transfer function is obtained through the Laplace transform, following the methodology of [19], which was developed for an undamped beam.
For identification of the model of flexible structures, frequency-domain methods are generally preferred over others in the time domain because the computational load required is lower, and for certain classes of systems, they are usually less sensitive to sensor noise and, consequently, may be more accurate [20]. Furthermore, identification algorithms in the frequency domain are traditionally founded on Levy’s method, which is a least-squares-based formulation whose results are not equally good at all frequencies, leading to the introduction of weights in some cases so as to handle this frequency dependence and improve the results [21]. Several methods can be found in the literature to identify a transfer function from a frequency response [22,23]. Some of them have been extended to fractional orders to describe some complicated real systems more adequately than integer-order models (see, e.g., [24] for further reading). In this respect, numerous studies have shown that many physical phenomena can be described more concisely and precisely by fractional-order models from frequency-domain data, such as in [25] for a permanent magnet synchronous motor, in [26] for the human arm dynamics, and in [27] for a fuel cell. Likewise, fractional-order models of flexible links were also identified in [28,29,30] in the frequency domain.
Given this motivation, the aim of this paper is threefold:
  • Measure, record, and characterize the intrinsic dynamics of the flexible link by extending the model proposed in [19] to fractional order, with the measured motor angle as model input and the measured motor coupling torque as output.
  • Consider different effects in the model in terms of damping. In particular, six different model structures will be considered, including internal and external damping, both of integer and fractional orders.
  • Analyze the variability in the identified models, as well as the performance of the fitting.
Then, the main contribution of this paper is to develop a general model, of fractional order, able to accurately describe the dynamics of flexible links including both internal damping, which is related to the link material, and external damping, which characterizes the environment in which the link is located. For illustration purposes, an application example is given for a carbon fiber antenna that moves freely in the air using frequency-domain identification methods.

3. Experimental Setup

The experimental setup is a version of one degree of freedom (DOF) of the robotic system with a single flexible link shown in Figure 3a, which was used as a sensing antenna in haptic applications (refer to, e.g., [31] for more details). It is a 2DOF platform, namely, it allows azimuth and elevation movement for the flexible link. From this system, 1DOF was removed (i.e., the elevation motor and its supporting structure) to allow only the link rotation on the horizontal plane. This modification aimed to avoid, as much as possible, the coupling of the resonant frequencies of the aforementioned motor structure set over the frequency response of the flexible link and, consequently, reduce noise in the experiments. Hence, the experimental platform used for this work is represented in Figure 3b, where it can be seen that a simpler component is used to attach the link to the azimuth motor.
To ensure stability, the structure has three legs made of stainless steel. The link is attached at one of its ends to one Harmonic Drive mini servo DC motor PMA-5A set, which includes zero-backlash reduction gears ( n = 100 ) and an incremental encoder. The encoder provides us with a precision of 7 × 10 5 rad outside the gear reduction. The system also has a force-torque (F-T) sensor, an ATI Mini 40, between the flexible link and the motor, allowing us to measure the coupling torque with a precision of 1.25 × 10 4 Nm. The structures that fix the F-T sensor to the motor and the link to the sensor are made of polylactide (PLA) with a 3D printer.
Table 1 shows the characteristics of the flexible link, which is made of carbon fiber. Here, L is the length of the link, ρ the linear mass density, E the Young’s modulus, and I the area moment of inertia about the bending axis.
Moreover, the motor is controlled in order to provide precise and fast motor positioning responses. The control structure proposed in [32] is used. It is a 2DOF proportional-integral-derivative (PID) controller with low-pass filters that ensures good trajectory tracking, compensates for disturbances such as unmodeled components of the friction, and is robust to parameter uncertainties.
For the data acquisition and control algorithms, the program LabVIEW 2010 is used. The sampling time is T s = 1 ms for data acquisition (measurements, control signals, and written data). In addition to this, MATLAB R2022b is utilized to implement the methods described below.

4. Dynamic Model

In this section, the behavior of a flexible sensing antenna is modeled using fractional-order models for the damping phenomena. With this aim, we obtain the exact transfer function between the motor angle (measured by the encoder) and the coupling torque of the antenna. This allows us to study only the dynamics of the antenna, without considering the motor behavior (which includes non-linear and time-varying friction components).
Therefore, we develop the dynamic model of a slewing flexible beam with fractional-order damping. Figure 4a details the scheme of the uniform flexible beam mentioned in the previous section. The following frames of reference are defined: an inertial frame ( X 0 , Y 0 ) and a non-inertial frame ( X , Y ) , which rotates attached to the motor. The rotation angle of ( X , Y ) with respect to ( X 0 , Y 0 ) is called the angle of the motor θ m ( t ) , w ( x , t ) is the transversal flexible deflection of the beam with respect to the frame ( X , Y ) , and Γ c o u p is the coupling torque between the beam and the motor.
The deflection w ( x , t ) is described by the Euler–Bernoulli beam theory [33] under the following hypotheses:
  • The deflection w ( x , t ) is small enough with respect to x to approximate arctan w ( x , t ) x by w ( x , t ) x .
  • The lateral deflection of the beam (out of the plane) is assumed to be zero.
  • The cross-sections of the slender beam remain plane and orthogonal to the beam axis even after deformation. Moreover, rotatory inertia and shear deformation can be neglected because our antenna is very slender.
Furthermore, damping phenomena are classically categorized as internal and external. Usually, when integer-order damping is taken into account, both dampings are incorporated into a single one. In other works, this single term is modeled as a fractional-order (see, e.g., [7]). In our model, both dampings are independently considered: internal damping is proposed to model the energy dissipated by the internal friction of the beam, which is related to the beam material, whereas external damping describes the beam vibrating in the environment.
Hence, the equation of motion is calculated from the equilibrium of moments in a differential segment of the beam. In this way, the two kinds of damping are included by using the appropriate force or moment [6]. Therefore, by considering the differential segment of the beam represented in Figure 4b, the following equation is obtained by the equilibrium of moments
M ( x , t ) x = Q ( x , t )
where M ( x , t ) is the bending moment, and Q ( x , t ) is the shear force [7].
The bending moment M ( x , t ) is correlated with the deflection w ( x , t ) . When internal damping is considered, as Kelvin–Voigt damping [6], the bending moment M ( x , t ) is the union of two moments that model the elastic and viscous behavior. Moreover, if a fractional-order model for viscoelastic materials [34] is used, the bending moment can be defined as
M ( x , t ) = E I + c d I α t α 2 w ( x , t ) x 2
letting c d and α R + ( 0 < α < 2 ) be the damping coefficient and the internal damping order, respectively.
The shear force Q ( x , t ) is calculated from the force equilibrium of the differential segment as
Q ( x , t ) x + q ( x , t ) = ρ 2 t 2 [ w ( x , t ) + x θ m ( t ) ]
where the right-hand term is the inertial force, and q ( x , t ) is the distributed force over the beam. Since we include external damping in the model, the distributed force q ( x , t ) is equal to the damped force exercised by the environment over the beam. This kind of damping is usually modeled as a damping force proportional to velocity when considering the air [6]. In our case, we assume it as a fractional-order damping force as follows
q ( x , t ) = c v λ t λ [ w ( x , t ) + x θ m ( t ) ]
where c v and λ R + ( 0 < λ < 2 ) are the damping coefficient and the external damping order, respectively. It should be highlighted that the total velocity of a differential segment of the beam is the addition of two terms: the deflection velocity and the linear velocity in x generated by the rotation of the motor.
Remark 1.
Typically, the fractional order α is confined to the range of 0 < α 1 , as noted in [11,34,35]. Nonetheless, we have elected to expand this definition to encompass the range of 0 < α < 2 , similarly to [7]. This assumption is also considered for the order λ.
Remark 2.
Fractional derivatives are defined in the Caputo sense [35].
Finally, the complete equation of motion, also called the Euler–Bernoulli equation, is obtained from Equations (1)–(4) as
E I w ( x , t ) + c d I α t α w ( x , t ) + ρ 2 t 2 w ( x , t ) + x θ m ( t ) + c v λ t λ w ( x , t ) + x θ m ( t ) = 0
where ( ) represents a spatial derivative with respect to x. In order to solve it, two initial conditions and four boundary conditions are needed. The four boundary conditions are the same as for a clamped-free beam, where the beam is fixed to the motor in x = 0 and it is free in x = L , i.e.:
w ( 0 , t ) = 0
w ( 0 , t ) = 0
w ( L , t ) = 0
w ( L , t ) = 0
Additionally, the relation between the coupling torque Γ c o u p ( t ) (which is equal to the moment at the link base) and the deflection w ( x , t ) is given by
Γ c o u p ( t ) = E I + c d I α t α w ( 0 , t )
Since we want to obtain the exact transfer function between Γ c o u p ( t ) and θ m ( t ) , the differential Equation (5), the boundary conditions (6)–(9), and Equation (10) are Laplace transformed with respect to the time derivative [19]. Then, Equation (5) results in
1 + ζ s α W ( x , s ) + ρ E I s 2 + γ s λ W ( x , s ) + x Θ m ( s ) = 0
where W ( x , s ) and Θ m ( s ) are the Laplace transforms of w ( x , t ) and θ m ( t ) , respectively, and the damping coefficients ζ = c d E and γ = c v E I are defined. Furthermore, making β 4 = ρ E I s 2 + γ s λ 1 + ζ s α , we have
W ( x , s ) β 4 W ( x , s ) + x Θ m ( s ) = 0
Additionally, the boundary conditions (6)–(9) become
W ( 0 , s ) = 0
W ( 0 , s ) = 0
W ( L , t ) = 0
W ( L , t ) = 0
and Equation (10) is
Γ c o u p ( s ) = E I 1 + ζ s α W ( 0 , s )
with Γ c o u p ( s ) the Laplace transform of Γ c o u p ( t ) .
In the first place, the system of partial differential equations formed by (11)–(15) is solved. This system can be expressed as a state-space representation with the space variable x as the independent one by
W ( x , s ) W ( x , s ) W ( x , s ) W ( x , s ) = 0 1 0 0 0 0 1 0 0 0 0 1 β 4 0 0 0 W ( x , s ) W ( x , s ) W ( x , s ) W ( x , s ) + 0 0 0 β 4 x Θ m ( s )
So, Equation (17) is
Z ( x , s ) = A ( s ) Z ( x , s ) + B ( s ) U ( x , s )
where we define
Z ( x , s ) = W ( x , s ) W ( x , s ) W ( x , s ) W ( x , s )
U ( x , t ) = x Θ m ( s )
A ( s ) = 0 1 0 0 0 0 1 0 0 0 0 1 β 4 0 0 0 , B ( s ) = 0 0 0 β 4
For the sake of simplifying the notation, the dependency of s on matrices A ( s ) and B ( s ) will be omitted from now on.
Despite Equation (18) not looking like a conventional state-space equation, it can be solved as a state-space model continuous and space-invariant. With that aim, we use the state-transition matrix. Then, the general solution is
Z ( x , s ) = Φ ( x , 0 ) Z ( 0 , s ) + 0 x Φ ( x , ϰ ) B U ( ϰ , s ) d ϰ
where the state-transition matrix, keeping in mind that A is constant for all x, is
Φ ( x , x 0 ) = e A ( x x 0 )
Moreover, the integral of Equation (19) can be solved taking into account that B is constant and using the form of U ( x , s ) . Then, the solution of Z ( x , s ) is
Z ( x , s ) = e A x Z ( 0 , s ) A 1 x I A 1 I e A x B Θ m ( s )
where the value of Z ( 0 , s ) needs to be known. With that purpose, the particular solution of (21) in x = L
Z ( L , s ) = e A L Z ( 0 , s ) A 1 L I A 1 I e A L B Θ m ( s )
is combined with the boundary conditions (12)–(15), resulting in a set of eight linear algebraic equations. By solving them, we can obtain the values of Z ( 0 , s ) and Z ( L , s ) .
Furthermore, the transition matrix can be evaluated using the inverse Laplace transform as Φ ( x ) = L 1 [ p I A 1 ] (converting a function dependent on p to another in the space domain). Since
p I A 1 = 1 p 4 β 4 p 3 p 2 p 1 β 4 p 3 p 2 p β 4 p β 4 p 3 p 2 β 4 p 2 β 4 p β 4 p 3 ,
the exponential matrix is
e A x = f ( x ) f ( x ) f ( x ) f ( x ) β 4 f ( x ) f ( x ) f ( x ) f ( x ) β 4 f ( x ) β 4 f ( x ) f ( x ) f ( x ) β 4 f ( x ) β 4 f ( x ) β 4 f ( x ) f ( x )
with  f ( x ) = 1 2 β 3 sinh ( β x ) sin ( β x ) .
Finally, taking Equation (16) with W ( 0 , s ) from the solution of  Z ( 0 , s ) , the exact transfer function between the coupling torque Γ c o u p ( s ) and the motor angle Θ m ( s ) can be obtained as:
G ( s ) = Γ c o u p ( s ) Θ m ( s ) = E I 1 + ζ s α β sin ( β L ) cosh ( β L ) cos ( β L ) sinh ( β L ) 1 + cos ( β L ) cosh ( β L )
The above transfer function is expressed in terms of transcendental functions of β , and differs from the result of [19] in the inclusion of the fractional damping, which affects the form of expressions (25) and β .

5. Methods

This section deals with the two-stage method used to identify the dynamics of the system. Firstly, details of the experiments are given, including the subsequent data processing carried out to estimate the frequency response of the flexible beam. Secondly, the methods used to identify the models of the flexible link are described. Finally, the indexes proposed to study the variability in the models are explained.

5.1. Experiment Description

Twenty experiments were carried out with our flexible link using LabVIEW. The reference applied to the motor control system, θ m * ( t ) , was a chirp signal with an amplitude of 0.025 rad, generated with frequencies distributed linearly in the [ 0.05 , 50 ]  Hz range during 120 s. The real motor angle, θ m ( t ) , and the coupling torque Γ c o u p ( t ) —which are, respectively, the system input and output according to (25)—were measured with a sampling frequency of 1 kHz. An example of signals involved in one experiment is illustrated in Figure 5.
To obtain the frequency response, an approach for estimating transfer functions of the form of (25) from spectral analysis was used. In particular, the measured frequency response of the link was computed by using Welch’s method—by means of the MATLAB function tfestimate—with 8 Hamming-windowed segments, 160,000 number of points, and 15% overlap in order to attain smooth frequency responses. It is important to remark that Welch’s method was chosen over the well-known fast Fourier transform (FFT) and others so as to gain a better estimation of the frequency response of the system in terms of noise. In particular, Welch’s method allows the division of the signals into segments, depending on the percentage of overlap defined, and to use windows, which improve frequency response estimation in the presence of noise. On the contrary, FFT obtains the spectral estimate over the entire signal, which produces poorer results in such a situation.

5.2. Identification Methods

Once the frequency response was found, G ( j ω p ) , p = 1 , , N , each fit was performed in MATLAB using a non-linear optimization procedure in the frequency domain by minimizing the following cost function
J = p = 1 N W t ( ω p ) G exp ( j ω p ) G ^ ( j ω p ) 2
where W t ( ω p ) was a normalized weighting function, G exp ( j ω p ) was the frequency response estimated from experiments, G ^ ( j ω p ) was the dynamic model frequency response calculated from (25) with the identified parameters, and N = 7907 was the total number of frequency response points resulting from the Welch’s method (i.e., frequency samples). The weighting function was defined as
W t ( ω p ) = 1 N in the case of constant weights 1 ω p i = 1 N 1 ω i 1 in the case of frequency - dependent weights
A detailed explanation of this function can be found in Appendix A. It should be remarked that other weighting functions were considered by giving more importance to resonances and anti-resonances of the frequency response, but the results were discarded as there were practically no differences with the case of constant weights.
Equation (26) was used as the cost function to identify six parameter sets corresponding to six different model structures of (25), depending on the type of internal and external damping considered, as follows:
  • Model with integer-order internal damping (henceforth referred to as type I).
  • Model with integer-order external damping (hereafter, type II).
  • Model with integer-order internal and external damping (type III).
  • Model with fractional-order internal damping (henceforth named type IV).
  • Model with fractional-order external damping (type V).
  • Model with fractional-order internal and external damping (type VI).
Note that model (25) has up to four parameters, namely, ζ , α , γ , and λ ; the remainder of the parameters, related to the flexible link properties, were invariant during identification and set to the values indicated in Section 3. The optimal values of the unknown parameters of each model were obtained to accurately fit the measured frequency responses by sweeping such parameters for different ranges and steps, which were gradually reduced by computation in order to find the global minimum of the cost function J: ζ [ 0.13 , 30.5 ] × 10 5 with a step of 0.01 × 10 5 , γ [ 0.012 , 0.45 ] with a step of 0.001 , whereas α and λ were varied within the interval ( 0 , 2 ) with a step of 0.001 . For clarity, Table 2 shows the values of the parameters that were fixed for each of the different models during identification. These models were obtained for the two forms of W t ( ω p ) ; consequently, twelve sets of parameters were actually identified for each measured frequency response.

5.3. Model Variability

The determination of the accuracy of the models plays a key role in system identification, in order to have a quality tag attached to them describing their uncertainty in terms of variability. Two types of variability were considered in this paper, related both to the model structure and to the nominal frequency response (obtained from the means of the fitted parameters ζ , α , γ , and λ ). The former refers to the statistical characteristics of the four identified parameters for the general case (i.e., ζ , α , γ , and λ ). The latter is measured through a group of performance indexes and refers to the goodness of the frequency response fit. The following statistics were computed: mean, median, standard deviation, minimum, and maximum. The performance indexes used were defined as follows:
  • Mean square error (MSE) per sampling frequency, defined as
    MSE = p = 1 N ( G p G ^ p ) 2 N
    where G p and G ^ p are the experimental and nominal frequency responses, respectively. It should be remarked that for W t ( ω p ) = 1 N in (26), this index coincides with the cost function J.
  • Weighted mean square error (WMSE), defined as
    WMSE = p = 1 N 1 ω p ( G p G ^ p ) 2 i = 1 N 1 ω i
    where ω p is the frequency in rad/s. This error allows more importance to be given to low frequencies. Similarly to MSE, this index coincides with the cost function J with frequency-dependent weights in (26).
  • Maximum deviation (MD), defined as
    MD = max | G p G ^ p |
  • Mean deviation of resonant peaks (MDR), which measures the mean deviation of each resonance mode of the frequency response of each experiment with respect to that of the nominal one in terms of both magnitude and frequency, defined as
    MDR ω q = k = 1 20 | ω q , k ω ^ q | 20 for q = 1 , 2 , 3
    MDR M q = k = 1 20 | | G | ω q , k | G ^ | ω q | 20 for q = 1 , 2 , 3
    where ω q , k and ω ^ q denote the frequency, in rad/s, of each resonant peak ( q = 1 , 2 , 3 ) for the experiments and the nominal model, respectively, whereas | G | ω q , k and | G ^ | ω q represent the magnitude, in dB, of the experimental and the nominal frequency responses, respectively, at the resonance frequency ω q .
It is important to highlight that the indexes MSE and MD are commonly used in the literature, whereas the other two were defined to study our system due to the importance of fitting the resonant peaks, especially the first one. In the habitual task of the sensing antenna, its behavior is mainly given by the low frequencies.

6. Results

This section contains the results obtained following the above-described methods. To find the model that accurately reflects each experimental data, the optimization method was applied separately to each of the twenty frequency responses experimentally measured (240 models in total, considering the six structures of the dynamic model (25) and the two forms of (26)). The mean values of the parameters of the resulting models for each structure are also referred to as “nominal”. Next, the variability in the results is explained in terms of both the model structure and the goodness of fit.

6.1. Variability in the Model Structure

Firstly, identification results for the six models considering the two forms of the weighting function are summarized in Table 3 (“const” refers to the case of W t ( ω p ) with constant weights, whereas “ ω -dep” represents frequency-dependent weights). The results are represented graphically in Figure 6, which compares the measured frequency responses (in gray color) with those obtained for the nominal models of each type considering the two forms of the weighting function W t ( ω p ) (the constant case in blue, and the frequency-dependent case in red). This comparison provides a clear understanding of the accuracy of the nominal models, and it can be achieved with a glance, allowing for efficient analysis. For the sake of conciseness, the actual parameters of all the 240 models are not reported. From these results, the following conclusions can be stated:
  • With respect to the measured responses plotted in Figure 6, although noisy at certain frequencies (mainly at very low and high frequencies and those corresponding to the anti-resonant peaks), we can confirm that the Welch’s method outperforms the FFT algorithm for estimating experimental frequency responses. Several zero-phase filters were tried on the time-domain experimental data and the Welch’s method configurations, but with poorer results than those shown. Likewise, taking into account the results shown in this figure, it can be said that qualitatively and looking at the fit of resonant peaks (in magnitude, phase, and frequency), the models whose frequency responses are most similar to the experimental ones are models type III and type VI (Figure 6c,f) (integer- and fractional-order models with two and four parameters, respectively). On the other hand, models type I and type IV (see Figure 6a,d) perform better at low frequencies and they excessively dampen high vibration modes because of parameter ζ (it should be recalled that this defines internal damping, which depends on the material of the flexible link). By contrast, models type II and type V (Figure 6b,e) are able to better adjust high frequencies due to the external damping parameter, γ . Obviously, these preliminary conclusions should be validated in a quantitative manner.
  • From Table 3, except for model type III, there are significant differences in the identified nominal model for each of the forms of weighting function W t ( ω p ) , even obtaining values for parameter ζ an order of magnitude higher when the frequency-dependent weights are considered. For the rest of the parameters, the mean values are reduced with the frequency-dependent weighting function, with parameter λ varying the least of these three, while α and γ are reduced by about 22% and 30%, respectively (except in the cases of the type V and VI models, where γ is practically unchanged or varies by about 75%, respectively). The differences can be better appreciated in Figure 6, excluding the case of model type III (Figure 6c), where the form of W t ( ω p ) in (26) has no influence (the frequency response of the models is the same). As can be observed, the main discrepancies are found in the resonance and anti-resonance peaks, both in magnitude and in phase (especially in the latter).
  • Considering the standard deviation values in Table 3, the parameter with the highest variability for all the identified models is ζ , with its value, except in the case of the model type III, being higher for the case of the frequency-dependent weighting function (note that this parameter does not appear in models type II and type V).
Figure 6. Experimental and identified frequency responses for models (a) type I, (b) type II, (c) type III, (d) type IV, (e) type V and (f) type VI. Here, (Fractalfract 07 00621 i003) is the frequency response of all the experiments, (Fractalfract 07 00621 i001) is the response of the nominal model obtained with constant weights, and (Fractalfract 07 00621 i002) is the one for the nominal model with frequency-dependent weights.
Figure 6. Experimental and identified frequency responses for models (a) type I, (b) type II, (c) type III, (d) type IV, (e) type V and (f) type VI. Here, (Fractalfract 07 00621 i003) is the frequency response of all the experiments, (Fractalfract 07 00621 i001) is the response of the nominal model obtained with constant weights, and (Fractalfract 07 00621 i002) is the one for the nominal model with frequency-dependent weights.
Fractalfract 07 00621 g006
Table 3. Model structure variability: statistics of the four identified model parameters. Note: “const” refers to the case of W t ( ω p ) with constant weights, whereas “ ω -dep” represents frequency-dependent weights.
Table 3. Model structure variability: statistics of the four identified model parameters. Note: “const” refers to the case of W t ( ω p ) with constant weights, whereas “ ω -dep” represents frequency-dependent weights.
Model α ζ ( × 10 4 ) λ γ ( × 10 1 )
const ω -depconst ω -depconst ω -depconst ω -dep
Type IMean--1.3310.75----
Mdn--1.32 10.71----
Std--0.090.34----
Min--1.1910.19----
Max--1.4611.42----
Type IIMean------1.971.39
Mdn------1.881.35
Std------0.290.20
Min------1.441.07
Max------2.581.82
Type IIIMean--0.540.54--0.200.20
Mdn--0.520.52--0.200.20
Std--0.080.06--0.040.01
Min--0.360.42--0.140.18
Max--0.700.65--0.290.23
Type IVMean0.8320.6792.4924.45----
Mdn0.8310.6782.4624.59----
Std0.0260.0090.231.01----
Min0.7680.6581.9722.45----
Max 0.876 0.697 2.95 26.68 ----
Type VMean---- 1.247 1.205 0.46 0.43
Mdn---- 1.252 1.209 0.39 0.43
Std---- 0.081 0.043 0.24 0.07
Min---- 1.035 1.110 0.20 0.31
Max---- 1.402 1.276 1.31 0.60
Type VIMean 1.362 0.996 0.09 0.71 0.616 0.594 2.25 0.55
Mdn 1.369 0.998 0.06 0.68 0.599 0.580 2.14 0.55
Std 0.147 0.057 0.14 0.20 0.082 0.028 1.04 0.04
Min 0.956 0.909 0.02 0.34 0.467 0.561 0.56 0.47
Max 1.605 1.124 0.64 1.13 0.805 0.658 4.43 0.61

6.2. Fitting Performance

Regarding the goodness of fit, two different fits can be evaluated by means of the indexes defined in Section 5: the first one is about the overall fit of the frequency response of the nominal models with respect to the experimental responses, whereas the second determines the performance of the nominal models with regard to the resonant peaks. Both are discussed next.

6.2.1. Fitting of Overall Frequency Response

Figure 7 shows the values of the MSE, WMSE, and MD for the 12 nominal models identified (from the six structures of model (25) and the two weighting functions) by means of box plots, to make the comparison easy. In the box plots, the green circles represent the mean values; red lines denote medians. The following issues can be concluded from these results:
  • Considering MSE (Figure 7a), the values obtained for all the models are lower for the constant weighting function, except for type III and type VI, for which there is practically no difference between the two forms of the function. Taking into account that this index is utilized as a cost function for identification (in this paper, the definition of MSE takes the form of the cost function (26) for W t ( ω p ) = 1 N ), constant weights are expected to show better results. However, it is remarked that there are no differences in model type III because the nominal models for the two weighting functions are the same. Likewise, when comparing all the identified models, models type III and type VI are the ones that perform better, with mean values of MSE of, respectively, 9.087 and 9.010 for constant weights, and 9.087 and 8.963 for frequency-dependent weights. The standard deviation is lower for the case of type VI. These results confirm the conclusion previously reached by looking at the frequency responses in graphical form. Again, it should be remarked that model type III is of integer order with two parameters, whereas model type VI is a fractional-order four-parameter dynamic model.
  • With respect to WMSE (Figure 7b), lower values for each model are obtained with frequency-dependent weights, as is expected (the definition of WMSE takes the form of the cost function (26) for frequency-dependent weights). Additionally, similar to the results obtained with MSE, the models with the worst performance are type I and type IV, whereas models type III and type VI have the best results, with mean values of 4.491 and 5.484 for constant weights, and 4.491 and 3.469 for frequency-dependent weights, respectively. In contrast with MSE, it is appreciated that frequency-dependent weights provide a better result than constant ones for models type I, IV, and VI. This fact agrees with the results shown in Figure 6f.
  • Where MD is concerned (Figure 7c), significant changes were found for the two weighting functions for models type I and IV, but for this index the worst results are obtained for W t ( ω p ) = 1 N .

6.2.2. Fitting of the Resonant Peaks

Since the accurate identification of the resonance modes of the sensing antenna is necessary, the performance of the nominal models on the resonant peaks is evaluated by employing the MDR indexes, which are reported in Table 4.
In general, there are significant differences between the models obtained with each of the two different weighting functions, both in frequency and magnitude. The first observation that can be drawn from the table is that, for the first resonant peak, and as expected, fitting results in terms of deviation are better with the frequency-dependent weights, especially in magnitude; in the case of frequency, there are no differences for models type I and type III. Secondly, analyzing each resonance individually for both frequency and magnitude, it can be said that the best results are obtained for:
  • First peak: in frequency, the best model is type VI with a mean deviation of 0 rad/s, followed by type II and type IV, with 0.039 rad/s. Additionally, model type I gives the best results in magnitude with a mean deviation of 0.288 dB. As expected, the best models, in this case, were all obtained with frequency-dependent weights, which by definition give more importance to the first resonant peak. Furthermore, except for type II, these models are dependent on internal damping that performs better at low frequencies.
  • Second peak: model type IV in frequency and model type VI in magnitude are the best, with constant and frequency-dependent weighting functions with a mean deviation of 0.112 rad/s and 0.956 dB, respectively.
  • Third peak: for frequency, model type VI with a constant weighting function (deviation value of 0.440 rad/s), and model type V with frequency-dependent weighting for magnitude (mean deviation of 0.865 dB) are the best fits. Hence, the models that best adjust the third peak are those that depend on fractional-order external damping.

6.3. Discussion about the Best Identified Models

Considering both the variability in the model structures and the performance of the fit, the model structures that best fit the measured frequency responses are the model with integer-order internal and external damping (type III) and the model with fractional-order internal and external damping (type VI), both with frequency-dependent weights. The following issues should be mentioned to sum up:
  • When comparing the two forms of cost function (26), the mean values of the damping parameters for model type III are practically equal for constant and frequency-dependent weights. However, the standard deviations of both ζ and γ are lower for the case of frequency-dependent weights (refer to Table 3). Furthermore, there is a slight improvement in the fit of the resonant peaks for frequency-dependent weights (see Table 4). Concerning model type VI, the standard deviations of the four damping parameters are also lower with frequency-dependent weights. Moreover, as was said, the difference between using the two weighting functions can be observed in Figure 6 and Figure 7.
  • For our system, using frequency-dependent weights for identification is preferable for the following reasons. Firstly, the behavior of the sensing antenna is mainly characterized by the first resonance mode. Hence, the accurate fit of the first resonant peak is preferred over the rest. Secondly, as is observed in experimental frequency responses in Figure 6, noise increases as frequency increases. While the noise at the anti-resonances has almost no influence, the noise produced in the third resonance peak substantially affects the fitting performance, as can be seen by comparing the MSE and WMSE indexes. MSE, which gives the same importance to all the frequencies, is more affected by this noise. Therefore, there are no significant differences between the identified models for the two weighting functions in terms of MSE, although their frequency responses are graphically quite different (see, e.g., Figure 6f).
  • Where variability in the model structure is concerned (Table 3), the standard deviations of parameters ζ and γ for model type VI (with percentages of 28% and 7% with respect to the mean value) are higher than those obtained for the parameters of type III (with 11% and 5%). Since model type VI has more parameters, it is expected to have more variability. Nevertheless, model type VI has a fractional-order internal damping that can be approximated to an integer order (mean and median of α are practically equal to 1). Meanwhile, the fractional-order external damping has a derivative order close to 0.5. If a model with fewer parameters than type VI is required, it can be set by a model with integer-order internal and fractional-order external dampings (with a number of three parameters). Furthermore, it could be said that integer order performs better in fitting internal damping than fractional order. Taking into account the considered experimental setup (carbon fiber beam in air), the results provide interesting information when compared with the literature (see, e.g., [34] and references therein): (a) α = 1 means that the beam can be described by classical Kelvin–Voigt damping; and (b) λ being close to 0.5 indicates that the interaction of the antenna with the environment can also be modeled as in other reference works.
  • Regarding the fitting of the overall frequency response (Figure 7), the difference between models type III and type VI in terms of MSE is slight: that of type VI is about 1% lower than the corresponding MSE of type III. However, this difference increases to 23% when considering WMSE. Since the noise of the third resonance peak affects MSE more than WMSE, it is better to use WMSE in our application.
  • Taking into account the deviation of resonant peaks in Table 4, the nominal model of type VI (obtained with the frequency-dependent weighting function) performs better than model type III in describing the first two resonance modes of the flexible link. The mean deviation of model type VI is 100% and 77% lower than type III in the frequency and magnitude of the first peak, respectively, and 14% and 44% for the second peak. However, the mean deviation of model type III is a 6% and 24% lower in frequency and magnitude, respectively, than the corresponding values of model type VI for the third peak.
In summary, although model type III has a better variability in its parameters, model type VI obtained with frequency-dependent weights is the one that best fits the experimental responses. Since damping for the antenna moving in the air has little influence on the time response, the interest in this study is to improve the accuracy of identification of its frequencies. The real frequencies of the system will depend on the physical characteristics (namely, link material and geometry) of the antenna and its damping. Therefore, the best model is a model with integer-order internal and fractional-order external dampings. As already mentioned, this result indicates that, for our antenna, the internal damping is described by the classical Kelvin–Voigt damping, while the external damping suggests a more complex behavior between the beam and the environment, which is sometimes modeled as non-linear, as in [5]. For a clearer comparison of the two models identified as the best, Table 5 contains the parameters of the corresponding nominal models, as well as their fitting performance (data extracted from Table 3, Figure 7b, and Table 4, respectively); their frequency responses are plotted together in Figure 8 (curves are from Figure 6).

7. Conclusions

This paper has addressed the dynamic modeling of a flexible antenna that is used for haptic navigation of mobile robots. The fast and vibration-free non-contact movement of the antenna is of utmost importance in order to obtain efficient haptic systems. This requires a high-performance controller, whose design must be based on a very precise dynamic model of the system.
We have developed a fractional-order general model that allows us to accurately describe the dynamics of flexible links involving internal and external damping, which are related to the link material and the environment where it is located, respectively. Likewise, we have investigated its adequacy for a carbon fiber antenna that moves freely in the air. In particular, six models with different combinations of dampings (of integer or fractional order) have been identified using frequency-domain methods from experiments conducted with our antenna. In total, 240 models have been obtained (20 experiments performed, six structures of the dynamic model, and two forms of the cost function considered for the non-linear optimization identification procedure). All the identified models have been assessed in terms of fitting accuracy to experimental curves and robustness (i.e., the variation in the estimated parameters when the identification is repeated several times).
After this process, we concluded that the best tradeoff between fitting accuracy and identification robustness is given by a model that includes an internal damping of integer order and an external damping of fractional order close to 0.5. Since damping has a small value in our antenna that moves freely in the air, the improvement attained by introducing fractional-order derivatives in the model has a low absolute value. However, the relative value of that improvement is noticeable.
Likewise, the proposed general model, of fractional order, will be valuable to accurately describe the dynamics of flexible links made of different materials that move within different environments by means of both the included internal and external damping, respectively. Our future work will be in such a direction: using the model to characterize the behavior of an antenna with different properties (e.g., made of carbon fiber and coated with a silicone pseudo-skin or of a polymer) immersed in water or other fluids.

Author Contributions

Conceptualization, M.I.H.-O. and V.F.-B.; methodology, M.I.H.-O., I.T., B.M.V. and V.F.-B.; software, M.I.H.-O. and I.T.; validation, M.I.H.-O. and I.T.; formal analysis, M.I.H.-O. and I.T.; investigation, M.I.H.-O. and I.T.; resources, V.F.-B. and B.M.V.; data curation, M.I.H.-O. and I.T.; writing—original draft preparation, M.I.H.-O. and I.T.; writing—review and editing, M.I.H.-O., V.F.-B. and B.M.V.; visualization, M.I.H.-O. and I.T.; supervision, V.F.-B. and B.M.V.; project administration, V.F.-B. and B.M.V.; funding acquisition, V.F.-B., I.T. and B.M.V. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Agencia Estatal de Investigación (Ministerio de Ciencia e Innovación) through the projects PID2019-111278RB-C21/AEI/10.13039/501100011033 and PID2019-111278RB-C22/AEI/10.13039/501100011033 and by the Agencia Estatal de Investigación and European Social Fund (ESF) through the scholarship PRE2020-095222 of the 2020 FPI Program.

Data Availability Statement

Not applicable.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

The following abbreviations are used in this manuscript:
DOFDegrees of freedom
FFTFast Fourier transform
maxMaximum value
MDMaximum deviation
mdnMedian
MDRMean deviation of resonant peaks
minMinimum value
MSEMean square error
PIDProportional-integral-derivative
stdStandard deviation
WMSEWeighted mean square error

Appendix A

Here, the two weighting functions of Equation (27) are detailed. The former is a constant function that sets the same importance for all the frequencies. In this case, Equation (26) is equal to the fitting performance MSE (Equation (28)). Usually, this cost function is taken to perform optimization. The latter is a frequency-dependent function, and its aim is to give more importance to the low frequencies of the system. Since high frequencies have less influence on the response of the flexible link, the first resonant peak should be a better fit.
To obtain the weighting functions, we start from the mean value theorem for definite integrals
ω 0 ω f e ( ω ) p ( ω ) d ω = e m ω 0 ω f p ( ω ) d ω
where ω 0 and ω f are, respectively, the initial and final frequencies, e ( ω ) is the defined error function, p ( ω ) is a generic weighting function, and e m is the mean value of the error function. The mean value e m is equivalent to the cost function J of Equation (26). So, Equation (A1) can be expressed as
e m = ω 0 ω f e ( ω ) p ^ ( ω ) d ω
by defining the normalized weighting function
p ^ ( ω ) = p ( ω ) ω 0 ω f p ( ω ) d ω
The equivalent equations for a discrete signal are
e m = i = 1 n e ( ω i ) p ^ ( ω i ) Δ ω
p ^ ( ω i ) = p ( ω i ) i = 1 n p ( ω i ) Δ ω
where n is the total number of points of the signal, and Δ ω is the constant interval of frequency variation.
Hence, Equation (A5) is equivalent to Equation (26) by setting W t ( ω i ) = p ^ ( ω i ) Δ ω . If p ( ω i ) = 1 , we obtain the constant weights of Equation (27).
Moreover, the frequency-dependent weights are defined to obtain the mean value of the error when a frequency logarithmic scale is considered. As the logarithmic scale is set when the frequency response is plotted, we seek to give the same weight to equally spaced frequencies in it. In order to obtain this weighting function, e ( ω ) is integrated with respect to another function g ( ω ) = log 10 ( ω ) , which is continuous and derivable by the Riemann–Stieltjes integral. The Riemann–Stieltjes integral is related to the standard Riemann integral by
ω 0 ω f e ( ω ) d g ( ω ) = ω 0 ω f e ( ω ) g ( ω ) d ω
where g ( ω ) = d g ( ω ) d ω . Therefore, the frequency-dependent weighting function can be defined with Equations (A1) and (A6) as p ( ω ) = g ( ω ) = 1 ω log ( 10 ) .

References

  1. Staudacher, E.; Gebhardt, M.; Durr, V. Antennal movements and mechanoreception: Neurobiology of active tactile sensors. Adv. Insect Physiol. 2005, 32, 49–205. [Google Scholar] [CrossRef]
  2. Prescott, T.J.; Pearson, M.J.; Mitchinson, B.; Sullivan, J.W.; Pipe, A.G. Whisking with robots from rat vibrissae to biomimetic technology for active touch. IEEE Robot. Autom. Mag. 2009, 16, 42–50. [Google Scholar] [CrossRef] [Green Version]
  3. Ueno, N.; Svinin, M.; Kaneko, M. Dynamic contact sensing by flexible beam. IEEE/ASME Trans. Mechatron. 1998, 3, 254–264. [Google Scholar] [CrossRef]
  4. Meirovitch, L. Analytical Methods in Vibrations; Macmillan Series in Applied Mechanics; Macmillan: New York, NY, USA, 1967. [Google Scholar]
  5. Baker, W.E.; Woolam, W.E.; Young, D. Air and internal damping of thin cantilever beams. Int. J. Mech. Sci. 1967, 9, 743–766. [Google Scholar] [CrossRef]
  6. Banks, H.T.; Inman, D.J. On Damping Mechanisms in Beams. J. Appl. Mech. 1991, 58, 716–723. [Google Scholar] [CrossRef]
  7. Li, M. Theory of Fractional Engineering Vibrations; De Gruyter: Berlin, Germany, 2021. [Google Scholar] [CrossRef]
  8. Mujumdar, A.; Tamhane, B.; Kurode, S. Fractional order modeling and control of a flexible manipulator using sliding modes. In Proceedings of the 2014 American Control Conference, Portland, OR, USA, 4–6 June 2014; pp. 2011–2016. [Google Scholar] [CrossRef]
  9. Singh, A.P.; Deb, D.; Agarwal, H. On selection of improved fractional model and control of different systems with experimental validation. Commun. Nonlinear Sci. Numer. Simul. 2019, 79, 104902. [Google Scholar] [CrossRef]
  10. Di Paola, M.; Heuer, R.; Pirrotta, A. Fractional visco-elastic Euler-Bernoulli beam. Int. J. Solids Struct. 2013, 50, 3505–3510. [Google Scholar] [CrossRef] [Green Version]
  11. Freundlich, J. Transient vibrations of a fractional Kelvin-Voigt viscoelastic cantilever beam with a tip mass and subjected to a base excitation. J. Sound Vib. 2019, 438, 99–115. [Google Scholar] [CrossRef]
  12. Zhang, G.; Wu, Z.; Li, Y. Nonlinear Dynamic Analysis of Fractional Damped Viscoelastic Beams. Int. J. Struct. Stab. Dyn. 2019, 19, 1950129. [Google Scholar] [CrossRef]
  13. Agrawal, O.P. Analytical Solution for Stochastic Response of a Fractionally Damped Beam. J. Vib. Acoust 2004, 126, 561–566. [Google Scholar] [CrossRef]
  14. Liang, Z.f.; Tang, X. Analytical solution of fractionally damped beam by Adomian decomposition method. Appl. Math. Mech. (Engl. Ed.) 2007, 28, 219–228. [Google Scholar] [CrossRef]
  15. Abu-Mallouh, R.; Abu-Alshaikh, I.; Zibdeh, H.S.; Ramadan, K. Response of Fractionally Damped Beams with General Boundary Conditions Subjected to Moving Loads. Shock Vib. 2012, 19, 321421. [Google Scholar] [CrossRef]
  16. Alkhaldi, H.; Abualshaikh, I.; Al-Rabadi, A. Vibration Control of Fractionally-Damped Beam Subjected to a Moving Vehicle and Attached to Fractionally-Damped Multiabsorbers. Adv. Math. Phys. 2013, 2013, 232160. [Google Scholar] [CrossRef] [Green Version]
  17. Behera, D.; Chakraverty, S. Numerical solution of fractionally damped beam by homotopy perturbation method. Open Phys. 2013, 11, 792–798. [Google Scholar] [CrossRef] [Green Version]
  18. Jena, R.; Chakraverty, S.; Jena, S. Dynamic response analysis of fractionally damped beams subjected to external loads using Homotopy Analysis Method (HAM). J. Appl. Comput. Mech. 2019, 5, 355–366. [Google Scholar] [CrossRef]
  19. Pota, H.R.; Alberts, T.E. Multivariable Transfer Functions for a Slewing Piezoelectric Laminate Beam. J. Dyn. Syst. Meas. Control 1995, 117, 352–359. [Google Scholar] [CrossRef]
  20. Tzes, A.P.; Yurkovich, S. Application and Comparison of On-Line Identification Methods for Flexible Manipulator Control. Int. J. Robot. Res. 1991, 10, 515–527. [Google Scholar] [CrossRef]
  21. Valério, D.; Ortigueira, M.D.; da Costa, J.S. Identifying a Transfer Function From a Frequency Response. J. Comput. Nonlinear Dyn. 2008, 3, 021207. [Google Scholar] [CrossRef]
  22. Pintelon, R.; Guillaume, P.; Rolain, Y.; Schoukens, J.; Van Hamme, H. Parametric identification of transfer functions in the frequency domain—A survey. IEEE Trans. Autom. Control 1994, 39, 2245–2260. [Google Scholar] [CrossRef] [Green Version]
  23. Ljung, L. System Identification: Theory for the User; Prentice-Hall, Inc.: Hoboken, NJ, USA, 1987. [Google Scholar]
  24. Valério, D.; da Costa, J.S. An Introduction to Fractional Control; The Institution of Engineering and Technology: Stevenage, UK, 2012; Chapter Fractional Identification; pp. 141–180. [Google Scholar] [CrossRef]
  25. Yu, W.; Luo, Y.; Chen, Y.; Pi, Y. Frequency domain modelling and control of fractional-order system for permanent magnet synchronous motor velocity servo system. IET Control Theory Appl. 2016, 10, 136–143. [Google Scholar] [CrossRef] [Green Version]
  26. Tejado, I.; Valério, D.; Pires, P.; Martins, J. Fractional order human arm dynamics with variability analyses. Mechatronics 2013, 23, 805–812. [Google Scholar] [CrossRef]
  27. Caponetto, R.; Matera, F.; Murgano, E.; Privitera, E.; Xibilia, M.G. Fuel Cell Fractional-Order Model via Electrochemical Impedance Spectroscopy. Fractal Fract. 2021, 5, 21. [Google Scholar] [CrossRef]
  28. Vinagre, B.M.; Feliu, V.; Feliu, J. Frequency domain identification of a flexible structure with piezoelectric actuators using irrational transfer function models. In Proceedings of the 37th IEEE Conference on Decision and Control (Cat. No.98CH36171), Tampa, FL, USA, 16–18 December 1998; Volume 2, pp. 1278–1280. [Google Scholar] [CrossRef]
  29. Feliu, V.; Feliu-Talegón, D.; Castillo-Berrio, C. Improved object detection using a robotic sensing antenna with vibration damping control. Sensors 2017, 17, 852. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  30. Traver, J.E.; Nuevo-Gallardo, C.; Rodríguez, P.; Tejado, I.; Vinagre, B.M. Modeling and Control of IPMC-Based Artificial Eukaryotic Flagellum Swimming Robot: Distributed Actuation. Algorithms 2022, 15, 181. [Google Scholar] [CrossRef]
  31. Mérida-Calvo, L.; Feliu-Talegón, D.; Vicente Feliu-Batlle, F. Improving the Detection of the Contact Point in Active Sensing Antennae by Processing Combined Static and Dynamic Information. Sensors 2021, 21, 1808. [Google Scholar] [CrossRef]
  32. Castillo-Berrio, C.F.; Feliu-Batlle, V. Vibration-free position control for a two degrees of freedom flexible-beam sensor. Mechatronics 2015, 27, 1–12. [Google Scholar] [CrossRef]
  33. Oñate, E. Slender Plane Beams. Euler-Bernoulli Theory. In Structural Analysis with the Finite Element Method Linear Statics: Volume 2. Beams, Plates and Shells; Springer: Dordrecht, The Netherlands, 2013. [Google Scholar] [CrossRef]
  34. Bagley, R.L.; Torvik, P.J. A Theoretical Basis for the Application of Fractional Calculus to Viscoelasticity. J. Rheol. 1983, 27, 201–210. [Google Scholar] [CrossRef]
  35. Mainardi, F. Fractional Calculus and Waves in Linear Viscoelasticity. An Introduction to Mathematical Models; Imperial College Press: London, UK, 2010. [Google Scholar]
Figure 1. Whiskers and antenna sensors in nature.
Figure 1. Whiskers and antenna sensors in nature.
Fractalfract 07 00621 g001
Figure 2. Prototype of mobile robot with sensing antenna system.
Figure 2. Prototype of mobile robot with sensing antenna system.
Fractalfract 07 00621 g002
Figure 3. Experimental platform: (a) 2DOF prototype; (b) 1DOF prototype (modified version for this work).
Figure 3. Experimental platform: (a) 2DOF prototype; (b) 1DOF prototype (modified version for this work).
Fractalfract 07 00621 g003
Figure 4. Dynamic model schemes: (a) slewing flexible beam scheme, (b) diagram of force in a differential segment of the link.
Figure 4. Dynamic model schemes: (a) slewing flexible beam scheme, (b) diagram of force in a differential segment of the link.
Fractalfract 07 00621 g004
Figure 5. Example of signals involved in one experiment: (a) motor response with (Fractalfract 07 00621 i004) reference angle and (Fractalfract 07 00621 i001) measured angle by the encoder; (b) measured coupling torque.
Figure 5. Example of signals involved in one experiment: (a) motor response with (Fractalfract 07 00621 i004) reference angle and (Fractalfract 07 00621 i001) measured angle by the encoder; (b) measured coupling torque.
Fractalfract 07 00621 g005
Figure 7. Boxplot of the fitting performance: (a) MSE, (b) WMSE, and (c) MD. Note: “const” refers to the case of W t ( ω p ) with constant weights, whereas “ ω -dep” represents time-dependent weights.
Figure 7. Boxplot of the fitting performance: (a) MSE, (b) WMSE, and (c) MD. Note: “const” refers to the case of W t ( ω p ) with constant weights, whereas “ ω -dep” represents time-dependent weights.
Fractalfract 07 00621 g007
Figure 8. Experimental and nominal frequency responses for models type III and type VI. Here, (Fractalfract 07 00621 i003) is the frequency response of all the experiments, (Fractalfract 07 00621 i001) is the response of the nominal model type III, and (Fractalfract 07 00621 i002) is the one for the nominal model type VI.
Figure 8. Experimental and nominal frequency responses for models type III and type VI. Here, (Fractalfract 07 00621 i003) is the frequency response of all the experiments, (Fractalfract 07 00621 i001) is the response of the nominal model type III, and (Fractalfract 07 00621 i002) is the one for the nominal model type VI.
Fractalfract 07 00621 g008
Table 1. Flexible link features.
Table 1. Flexible link features.
ParameterValueUnit
L 0.98 m
ρ 4.7 × 10 3 kg/m
E 1.15 × 10 11 N/m2
I 7.85 × 10 13 m4
Table 2. Values of parameters ζ , α , γ , and λ that were fixed for the six considered models for the identification procedure.
Table 2. Values of parameters ζ , α , γ , and λ that were fixed for the six considered models for the identification procedure.
ParameterModel
Type IType IIType IIIType IVType VType VI
α 1010
ζ 00
λ 0110
γ 00
Table 4. Mean deviation in resonant peaks.
Table 4. Mean deviation in resonant peaks.
Resonant PeakWeighting FunctionModel
Type IType IIType IIIType IVType VType VI
Frequency (rad/s)
Firstconst 0.118 0.157 0.118 0.118 0.471 0.432
ω -dep 0.118 0.039 0.118 0.039 0.353 0.000
Secondconst 0.269 0.151 0.269 0.112 0.858 0.190
ω -dep 0.673 0.230 0.269 1.459 0.662 0.230
Thirdconst 2.250 1.818 1.818 3.114 0.997 0.440
ω -dep 20.471 1.779 1.818 13.717 1.268 1.936
Magnitude (dB)
Firstconst 17.488 22.922 3.635 15.788 15.778 13.399
ω -dep 0.288 20.069 3.571 1.809 14.473 0.812
Secondconst 4.883 14.102 1.746 3.136 10.822 2.642
ω -dep 22.800 11.102 1.704 15.492 8.787 0.956
Thirdconst 12.245 3.670 5.195 8.969 2.613 6.345
ω -dep 27.427 0.993 5.167 19.262 0.865 6.790
Table 5. Parameters and performance of the nominal models type III and type VI with frequency-dependent weights (best models).
Table 5. Parameters and performance of the nominal models type III and type VI with frequency-dependent weights (best models).
ModelParameterPerformance Index
α ζ ( × 10 5 ) λ γ ( × 10 2 )WMSE MDR ω 1 (rad/s) MDR M 1 (dB) MDR ω 2 (rad/s) MDR M 2 (dB) MDR ω 3 (rad/s) MDR M 3 (dB)
Type III1 5.4 1 2.0 4.491 0.118 3.571 0.269 1.704 1.818 5.167
Type VI 0.996 7.1 0.594 5.5 3.469 0.000 0.812 0.230 0.956 1.936 6.790
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Haro-Olmo, M.I.; Tejado, I.; Vinagre, B.M.; Feliu-Batlle, V. Fractional-Order Models of Damping Phenomena in a Flexible Sensing Antenna Used for Haptic Robot Navigation. Fractal Fract. 2023, 7, 621. https://0-doi-org.brum.beds.ac.uk/10.3390/fractalfract7080621

AMA Style

Haro-Olmo MI, Tejado I, Vinagre BM, Feliu-Batlle V. Fractional-Order Models of Damping Phenomena in a Flexible Sensing Antenna Used for Haptic Robot Navigation. Fractal and Fractional. 2023; 7(8):621. https://0-doi-org.brum.beds.ac.uk/10.3390/fractalfract7080621

Chicago/Turabian Style

Haro-Olmo, María Isabel, Inés Tejado, Blas M. Vinagre, and Vicente Feliu-Batlle. 2023. "Fractional-Order Models of Damping Phenomena in a Flexible Sensing Antenna Used for Haptic Robot Navigation" Fractal and Fractional 7, no. 8: 621. https://0-doi-org.brum.beds.ac.uk/10.3390/fractalfract7080621

Article Metrics

Back to TopTop