Measuring the maximal Lyapunov exponent from data

One of the main goals of the dynamic systems analysis is to obtain the maximal Lyapunov exponent from a time series generated by the system. For this purpose, several algorithms have been implemented: Wolf et.al (1985), Sano and Sawada (1986) and Rosenstein el.al.

(1993). These methods have in common the main numerical procedure for the calculation of the maximal Lyapunov exponent from a time series. A reconstruction of the phase space and the test of the exponential divergence of nearby trajectories by means of an average

5The bounded orbits corresponding to aC¹ class differential equation over the plane are periodic orbits, or (asymptotically) tend to periodic orbits or end in a fixed point [75].

value of those divergence along the whole data set.

Given a time series with values{x₁, x2, . . . , xn}, the trajectory is expressed in the phase space reconstructed by delays as a matrix X with dimension M×m :

X=

withmbeing the optimal dimension for the phase space reconstructed by delays or embed-ding space, M the number of initial points that belong to the first delay interval (as shown in figure 2-10) andν the temporary delay nondimensionalized with the sampling period of the original time series (ifJ is the value of the time delay, andf_s is the sampling frequency of the time series, the relation ν = J ·fs is verified). In order to build the matrix, it is necessary to consider a duration of the time series as a number N of points such as N is multiple of the optimal dimension of the reconstructed phase space. In this way the relation N =M×m is verified.

Figure 2-10: Representation of the time series points and their correspondence with the elements of the matrix M, equation (2.75).

Only with the time series, the whole state of the system is defined: e.g. in the i−th discrete time, the system state is given by

X_i =

x_i, x_i+ν, . . . , x_i+(m−1)ν

(2.76) Notice that every vector-row belonging to the trajectory matrix is a point in the recon-structed phase space. The numerical method for the maximal Lyapunov exponent calcu-lation consists in taking one of these points X_n0 as a reference, selecting all the points in a sphere of radius centered in Xn0, and computing the average of the distances between

the point X_n0 and all these nearby points as a function of the relative time (expressed as a delay multiple). The logarithm of the averaged value of the ratio between the distance from the reference point to each point in the sphere and that distance in the initial time is a measure that contains all the deterministic fluctuations due to either the projection over the reconstructed phase space or the inherent dynamics of the system. Finally, an average of those logarithmic values must be done by taking as reference Xn0 the other points (an arbitrary set of N points) of the time series. The values of last average depend on the discrete time ∆nof the time series:

λ₁(∆n) = 1

where U(Xn0) is the set of points Xn belonging to the the sphere of radius centered in the reference point X_n0, d_Xn(0) is the euclidean distance in the initial time between the reference point Xn0 and the generic point Xn belonging toU(Xn0), and dXn(∆n) is the distance between these two generic points after a time lapse ∆n (note that this time lapse implies a displacement of number of points ∆n in the reconstructed-by-delays matrix of the system). λ₁ represents an estimation of the maximal Lyapunov exponent. From this basic procedure each algorithm introduces slight variations to improve the convergence of the numerical method, and to make it more robust and reliable.

A priori, the optimal value of the radius is unknown, and the problem of the presence of noise in the time series data must be considered. If the noise level is higher than the radiussome points may be falsely considered to belong to the topological sphere. On the other hand, if the considered dimension m for the reconstructed phase space is big enough to avoid the intersection of trajectories, small variations ofm wont affect the result.

An improved version of the Wolf algorithm implemented in Matlab is used in this thesis to perform the estimations of the maximal Lyapunov exponent, as a way to explore the presence of chaotic behavior in the generated time series. The main feature of the Wolf algorithm is the that a single neighbor is followed and repeatedly replaced when its sepa-ration from the reference trajectory grows beyond a certain limit [15], [77]. To corroborate its proper performance, several time series, obtained from the numerical resolution of the Duffing equation, have been used to test the algorithm, and its results have been compared with those provided by the Govorukhin algorithm. The Govorukhin algorithm, developed

by V. Govorukhin is a free source code which provides the Lyapunov spectrum from the ODE -system instead of the time-series. For each set of parameters corresponding to the equation parameters, the Govorukhin algorithm obtains the evolution of the three Lyapunov exponents of the Duffing system. The good agreement between both methods, the improved and Matlab-implemented Wolf algorithm and the Govorukhin algorithm (shown in figures 2-11 and 2-12, for either linear and chaotic Duffing systems), proves the accuracy of the first one, which can be used for time-series obtained from any source. The comparison is per-formed between the maximal Lyapunov exponent provided by the Govorukhin method and the Lyapunov value (corresponding to the maximal exponent) given by the Wolf algorithm.

Figure 2-11: Lyapunov exponent of a Duffing linear system obtained from a) its time series with Wolf algorithm and b) its ODE system with Govorukhin algorithm.

Figure 2-12: Lyapunov exponent of a Duffing chaotic system obtained from a) its time series with Wolf algorithm and b) its ODE system with Govorukhin algorithm.

Chapter 3

MEMS resonators: theory and fabrication

In this chapter, the operating principles and the fundamental parameters of electrostatically actuated cc-beam resonators with capacitive readout are described. The used analytical approach is the one degree of freedom (1DOF), obtained from the Galerkin discretization and order-reduction method. Taking this approach as starting point, an equation that govern the dynamical behavior of a generic cc-beam is obtained and all its terms are modeled with analytical expressions. In the same way, the electrostatic actuation and the capacitive readout principles are exposed and their influence over the MEMS performance is analyzed.

The second part of this chapter relates the used fabrication and integration techniques and approaches in the framework of CMOS technologies, considered in this work.

3.1 MEMS modeling

The electromechanical system under consideration consists in a clamped-clamped beam placed symmetrically and in plane between two electrodes in which electrostatic actuation and capacitive readout may be respectively performed. Figure 3-1 depicts the scheme of the system and shows its geometrical parameters.