On the Computation of Floquet Multipliers for Periodic Solution in Piecewise-smooth Dynamical System

The floquet multiplier is one of the most important indicators for the stability and bifurcation analysis for periodic solutions in nonlinear dynamical systems. Different from the well-established Floquet theory for the perturbation systems of smooth systems, much less has been understood in its counterpart for non-smooth systems. Here in this paper, we will report an unusual and interesting feature of the Floquet multipliers for piecewise-smooth dynamical systems. When the initial condition of the periodic solution is located at the boundary splitting the solution domain, the multipliers would be calculated falsely in certain circumstances, respectively, by a saltation matrix method or a direct numerical integration for the perturbation system. We elucidate the origin of the fake multipliers through perturbation analysis, and furthermore suggest an effective manner to avoid the miscalculation. This finding would be of fundamental significance to both the real-world applications and theory establishment of the Floquet theory in non-smooth systems


Introduction
As one kind of the most common responses arising in nonlinear dynamical systems, periodic solution has always been one of the centers of research in controlling chaotic motions [1][2] and stabilizing/destabilizing certain system responses [3][4][5][6].For a long period of time, the Floquet multiplier theory has been an indispensable tool in the stability and bifurcation analysis for periodic responses [7][8][9][10][11].This is an important and fundamental issue, as the stability evaluation is not only associated with system analysis [7], but also one of the cornerstones for realizing the stabilizing [4] or destabilizing [6] of certain dynamic behaviors.
The classical Floquet theory, however, was developed based on the assumption of system smoothness.For non-smooth systems, peculiar phenomena can possibly appear such as the jump of multipliers [12].Though important contributions have been made to better understand or extend this theory to non-smooth systems for more than a half of century [13][14][15], to the best of our knowledge there is yet a lack of rigorous analysis for the complete theoretical framework.
In this paper we will present an interesting phenomenon that, under certain conditions the multipliers would be obtained falsely for non-smooth dynamical systems, which are called as fake multipliers for convenience.To introduce the Floquet theory, first we consider a smooth autonomous system ) ( where the superscript denotes the differentiation with respect to time t.The n-dimensional vector, ) (x f , is differentiable to the first-order at least. .To check the stability of the periodic solution, according to the Floquet theory, is to evaluate the orbit convergence or divergence of a neighboring solution from the periodic orbit.Substitute an initial perturbation into (1) with , we obtain the variation equation [16] x with j e as a unit vector with the j-th element as 1, we can get the state at . With the attained solutions, we construct the monodromy matrix as )] ( , ), ( [ , whose eigenvalues are called as the Floquet multipliers.
For the sake of conciseness, we take an illustrative example [18] with an additional cubic nonlinearity [ ] , the coefficient vectors v and 3 v , the coefficient matrix A , and a given parameter μ .Herein, the superscript " T " denotes the transpose of a vector or matrix and ) sgn( ⋅ denotes the sign function.We firt consider a piecewise linear system (i.e., 0 v =

3
) and rewrite Eq. (4) as here the boundary line separates the phase plane as the upper ( ) and the lower ( . And, a more complicated example with a cubic nonlinearity will be considered at the end of this article. As mentioned above, the Floquet theory was established based on the system smoothness.As the considered system is non-smooth, the variation equation cannot be deduced as the right hand side is indifferentiable.Alternatively, the variation equation can be obtained at each region.As it is shown in Fig. 1(a), without loss of generality we assume the initial state ( ini x ) of the periodic solution is located at the lower region with is the second entry of ini x .The periodic curve crosses the boundary at L x for the first time after 1 t and at R x for the second time after 2 t , then it goes back to ini x after T .With these denotations, we obtain the variation equation at different regions Integrating Eq. ( 6) over ) , 0 ( 1 t under the same initial conditions utilized for (2), an matrix can be obtained as )] ( ), ( [ . Due to the jump of the fundamental solution matrix in discontinuous systems, Aizerman and Gantmakher [13] constructed an additional matrix, named as saltation matrix, to compensate this discontinuity.More details about the construction of the saltation matrix can be found in [16].For the considered system, specifically the solution matrix of ( 6) should be updated as with L S as the saltation matrix associated with the left intersection point, that is with a scalar and a row vector [16].Take 1 L M S as the updated initial conditions and integrate (6) over ) , ( 21 t t , we obtain . Finally, we obtain the monodromy matrix at .During the integrating process of the variation equation, the solution matrix has been premultiplied twice by the corresponding saltation matrix.Surprisingly, the multipliers are obtained to be 1.2177 and 0.0057 if the premultiplication is done only once when shown in Fig. 1(b), and to be 1.1853 and 0.0011 when L ini x x = .This indicates the appearance of fake Floquet multipliers.It should be pointed out that there is no multiplier equal to 1, which violates the Floquet theory that there is a unit multiplier for periodic solutions of autonomous systems [28].In practice, the leading multiplier larger than 1 will lead to a misjudgement on the stability of the solution.

Direct Numerical Integration for Multipliers
In this section we will present a direct numerical integration algorithm to approximate the multipliers, in order to study intensively the appearance of the fake multipliers.Recall that the original idea of the multipliers is to evaluate the departure or approaching of a perturbative solution to the periodic orbit.Assume the periodic solution y and 2 y will end with the returning states after one period, denoted by 1 z and 2 z respectively. .As ini M and ret M are dependent upon ini x , 1 y and 2 y , we calculate the monodromy matrix asymptotically as the limit of Since there is no analytical technique to solve the above limit, we approximate it by selecting the initial perturbations close enough to the initial state.When the initial state is .The numerical results converge to the true multipliers provided above by the method based on the saltation matrix, as long as the perturbative curves are close enough to the periodic orbit.To highlight the main results from the direct calculation of the multipliers, a more thorough analysis will be made at the end of this section.
Using the asymptotic approach, we come to the following conclusions about the influences of the initial state and the perturbations on the eigenvalues of the asymptotic matrix (1) The multipliers can be calculated correctly as long as the initial state is not located at the boundary.( 2) If the initial state is at the boundary, only when the initial perturbations and their returning states are all located at the same region, can the true multipliers be obtained.to indicate which regions i y and i z are located at, for the purpose of discussing the influences of the initial perturbations and the returning states.The sign "+" means the corresponding state is located at the upper region, and " − " for the lower region.For example, ) , ( + + represents that both i y and i z are located at the upper region; and ) , ( − + indicates that i y and i z are located at the upper and lower region, respectively.And the rest is defined in the same manner.The second point raised above is illustrated by Fig. 2(a), indicating that the true multipliers are obtained as long as the initial perturbations and their returning states are located at the same region.
From Fig. 2(b) we observe that the eigenvalues of ) y y J(x converge to constants other than the true multipliers, when the initial perturbations are at one region and their returning states at the other.Interestingly, the fake multipliers provided in the case of ) , ( + − coincide with the results attained by the method based on the saltation matrix.Aside from the cases that the initial perturbations are at the same region and so are the returning states, the eigenvalues do not converge to any constants, as shown in Fig. 2(c).
Due to the discontinuity in the right hand side of system (5), there are two different tangential vectors at the initial state, which can be normalized to be at the lower region, and at the upper region.For the case of ) , ( + + , the eigenvector corresponding to the unit multiplier is parallel to 2 γ , and the counterpart for ) , ( − − is parallel to 1 γ .Whereas for the fake multipliers, the eigenvector corresponding to the leading multiplier in the case of ) , ( − + is parallel to 1 γ , and to 2 γ in the case of ) , ( + − .In a word, the eigenvector for the leading multiplier (either true or fake) is parallel to the tangential vector at the region where the returning states are located at.
To examine the convergence of the multipliers attained by the direct calculation, we recall that and denote i i q x y ε + = ini with i q as scaling vectors and 1 << ε as a scaling factor.By this manner the initial perturbation matrix is given as [ ]  (10) Notice that, both the asymptotic monodromy matrix + + above method we obtain the multipliers as 1.0000 and 0.0012 with the initial state be randomly selected as T line, we present in Fig.1(b) the periodic orbit as well as two perturbative curves.The curves starting from initial states 1

Figure 1 .
Figure 1.(a): Phase planes of the periodic solution and perturbative curves for system (5) with 2 . 1 = μ lower region, the multipliers are obtained to be 0

Figure 3 .
Figure 3. Convergence of the leading multiplier provided by the monodromy matrix (7) versus the scaling factor ε at T R ] 0.91454492 [ μ − = x, with the scaling vectors as