A Time-Domain Asymptotic Approach to Predict Saddle-Node and Period Doubling Bifurcations in Pulse Width Modulated Piecewise Linear Systems

In this paper closed-form conditions for predicting the boundary of period-doubling (PD) bifurcation or saddle-node (SN) bifurcation in a class of PWM piecewise linear systems are obtained from a time-domain asymptotic approach. Examples of switched system considered in this study are switching dc-dc power electronics converters, temperature control systems and hydraulic valve control systems among others. These conditions are obtained from the steady-state discrete-time model using an asymptotic approach without resorting to frequency-domain Fourier analysis and without using the monodromy or the Jacobian matrix of the discrete-time model as it was recently reported in the existing literature on this topic. The availability of such design-oriented boundary expressions allows to understand the effect of the different parameters of the system upon its stability and its dynamical behavior.


Introduction
Switched systems constitute a special class of nonlinear systems [1] and arise often in many practical engineering systems when some switching elements such as switches or diodes, block with dead-zone, saturated amplifiers, relays and comparators in electrical systems are present. This is also the case of mechanical systems where impacts or nonsmooth friction take place. A particular class of switched systems are those characterized by linear differential equations between switching events. These systems are called therefore piecewise linear (PWL) or piecewise affine (PWA) systems [2]. Most of the PWL systems studied in the literature are characterized by switching among linear subsystems when certain time-varying and T − periodic boundaries in the state space are reached. This is the case of Pulse Width Modulation (PWM) systems like switching dc-dc power converters [3], [4], [5], [6], [9], [10], dc-ac inverters [11], temperature control systems [12], switched capacitor networks and chaos generators [13] and hydraulic and fluid valve drivers [14], [15]. In steady-state, during a switching period of length T , a trajectory of these systems starts at time instant nT and is described by the vector field f 1 (x) = A 1 x + B 1 u, intersects a switching boundary described by the equation σ(x(t), t) := Fx(t) − r(t) = 0 at switching instant t s = DT , and then goes to another linear system described by the vector field f 2 (x) = A 2 x + B 2 u, where r is a time-varying T −periodic external signal, x ∈ R n is the vector of the state variables, n is the order of the system A i ∈ R n×n and B i ∈ R n×m , i = 1, 2 are the system state matrices for phase i (i = 1, 2) and u ∈ R m is the vector of the system inputs in both the plant and controller, m being the number of the external inputs to the system which are supposed to be constant within a switching cycle. The ratio D = t s /T is called the duty cycle. The system can be modeled as follows a e-mail: abdelali.elaroudi@urv.catẋ For each phase, the system equations are linear and time-invariant and can be solved in closed-form. The trajectory x(t) at time t of the system starting for an initial condition x(t 0 ) at time instant t 0 can be expressed as follows Nonlinearity arises from the feedback which imposes the constraint σ(x(DT ), DT ) := Fx(DT ) − r(DT ) = 0 relating the duty cycle D nonlinearly and in general implicitly to the vector of state variables x of the system. Despite their engineering use, one of the main drawbacks of PWL and PWM systems is this nonlinearity making them prone to exhibit a large variety of complex dynamic and undesired behaviors [5], [12], [13]. Although each subsystem is linear and its describing differential equations can be solved in closed-form, the dynamics of the complete switched system is highly nonlinear and its analysis requires sophisticated computational tools [16]. The dynamical behavior and stability analysis of this kind of systems can either be tackled by long-time integration of the continuous-time switched model, discrete-time model and its Jacobian matrix or Floquet theory with Fillipov technique to compute the mondromy matrix [17], [18], [19]. After obtaining the Jacobian or monodromy matrix, critical boundary conditions for some singularities like saddlenode (SN) bifurcation or period-doubling (PD) can be obtained by imposing in the characteristic equation in the zdomain that one of the eigenvalues is equal to +1 or −1 respectively [4]. Another approach recently used in [6] for locating these boundaries which was wrongly called harmonic balance is by expanding the feedback signal into a Fourier series to obtain the steady-state trajectory in certain periodic regimes and imposing critical conditions for the occurrence of the corresponding singularities like PD and SN bifurcations 1 . With that approach an effort to transform the results from the Fourier frequency-domain into the time-domain has still to be done. In [6], the transformation from the Fourier frequency-domain to the timedomain is based on elementary partial fraction decomposition after defining some elementary cases of the system transfer function in the s−domain and listing them in the form of tables. However, this transfer function cannot be directly defined for systems with A 2 A 1 making the approach only applicable for a limited class of PWM systems like the ones considered in [6]. In particular those that can be be formulated in the form of a linear system and a square-wave siganl generated by a comparator like the PWM process. In [6] the author applied the approach based on the Fourier series expansion to a system with A 2 A 1 only by making some approximations leading finally to A 1 = A 2 . Another different type of approximation leading to the same consequence A 1 = A 2 has been used and justified in [9] for this kind of systems.
In [20], the Poisson sum formulae and some related Fourier series properties have been used to transform the condition for PD occurrence derived in [6] from the Fourier frequency-domain to a matrix-form state-space time-domain condition. For design purpose, these matrix-form expressions can be expressed in scalar form by using some practical assumptions [3]. We will see in this paper that a steadystate analysis of the trajectory in the time-domain will lead to the same results without need to use the Fourier series expansion and without having to perform any transformation nor needing the calculation of the Jacobian or the monodromy matrix.
The rest of the paper is organized as follows. Section 2 presents the derivation of a closed-form expression for predicting SN singularity in a feedback PWM system with PWL nonlinearity starting from the expression of the steadystate T −periodic orbit in the time-domain at the switching instant and imposing the constraint σ(x(DT ), DT ) = 0 due to the feedback. Subsequently, the work derives in Section 3 a similar condition for predicting PD bifurcation in this kind of systems by imposing a boundary condition in the time-domain between T -periodic and 2T −periodic orbits. Finally some concluding remarks are drawn in the last section.

Asymptotic methods for bifurcation boundary prediction
Let x ss (0) be the steady-state value of the periodic orbit of the system at the beginning of the period and x ss (DT ) be the steady-state value of this orbit at time instant DT . The 1 The author of [6] call this type of analysis a Harmonic Balance approach but one cannot see where this Harmonic Balance was performed. One can observe that equating the steady-state feedback signal Fx expressed in Fourier series to the external periodic signal at the switching instant cannot be called a Harmonic Balance as defined by experts in this field (See [7], [8] for example). Note that the Harmonic Balance approach implies substituting the Fourier series expression in the dynamical model of the original system and equating terms that have the same frequency.

r(t)
Waveforms of the T −periodic external signal r(t) and the control signal Fx ss (t) at T −periodic regime in steady-state.
system is forced periodically and synchronously to the first phase (i = 1) characterized by the vector field f 1 while it is switched to the second phase (i = 2) characterized by the vector field f 2 whenever Fx(t) − r(t) = 0 (See Fig. 1).
Therefore, in steady-state, according to (3), the vector of state variables at the beginning of the switching period is given by ( Fig.  1) In turn, the vector of state variables at the switching time DT within the switching period can be expressed as follows Putting (4) in (5), one obtains For more simplicity, let time t be normalized with respect to the switching period T in such a way that the switching duration t s coincides with the duty cycle D. Therefore the state variables x ss (D) at the instant D in steady-state is given by Let Φ = Φ 1 Φ 2 and Ψ = Φ 1 Ψ 2 + Ψ 1 . Therefore (7) can be simplified as follows where the matrix (I−Φ) is assumed to be nonsingular. If an integrator exists in the system, this matrix can be singular but some tricks exist to avoid this singularity.

Derivation of saddle-node (SN) bifurcation boundary in the time-domain
As mentioned before, the feedback loop imposes the following constraint between the state variables x(D) and the duty cycle D For the sake of simplicity of notation let x ss (D) = x(D). SN bifurcation or tangent bifurcation is a local bifurcation in which two solutions of a continuous system collide and annihilate each other. In discrete dynamical systems, the same bifurcation is often instead called a fold bifurcation or blue skies bifurcation in reference to the sudden creation of two fixed points. When the discrete-time system is obtained by sampling a continuous-time system, these fixed points correspond to periodic orbits of the original system and the phenomenon can also be called cyclic fold bifurcation. Fig. 2 illustrates the occurrence of a SN bifurcation in a nonlinear system.

Bifurcation parameter
T -periodic node x T -periodic saddle In PWM systems considered in this study, at the boundary of a SN bifurcation, there is a tangency between Fx(t) and r(t) in such a way that two solutions of (9) coalesce and disappear. Therefore, from (9), the following equality holds at this critical point Let σ e (D) = ∂r/∂D be the slope of the external T −periodic signal r(t) at time instant D. Therefore (10) becomes The derivative of the left side of (11) can be obtained by using (6) (note that x(D) = x ss (DT )) and differentiating the involved matrix functions. Let us calculate the derivative ∂x(D)/∂D. Using (6) one obtains Using known chain rules for matrix derivative, (12) can be written as follows Putting the term (I − Φ) −1 as a common factor results in the following equation Then by using (8), (14) becomes The derivative of the involved matrix function ∂Φ(D)/∂D can be calculated as follows After some algebra, the derivative ∂Ψ(D)/∂D can be expressed in the following form Let ΔA = A 1 − A 2 and ΔB = B 1 − B 2 . Substituting (16) and (17) in (15), the critical boundary condition for SN bifurcation boundary in (11) becomes where σ e,SN (D), the critical slope of the external function r(t) for SN bifurcation occurrence, can be expressed by Using (4), (19) becomes Taking into account that , the critical value of the slope of the external T −periodic function at the boundary of a SN bifurcation is where . It has to be mentioned here that in [4], [10] a slightly differently expressed condition has been obtained for the same boundary condition which is reported and adapted here for comparison Although apparently the condition (22) derived in [4] and [10] and that in (21) do not coincide, they just happen to be the same conditions but expressed differently. Note however that the expression (21) is simpler than the one derived in [4] and [10]. It must be pointed out that in [6], the steady-state asymptotic approach was only applied to a dc-dc switching buck converter for which A 1 = A 2 obtaining the stability condition only for similar particular cases. For obtaining the same equivalent stability boundaries for systems with A 1 A 2 , two different approaches have been used in [10], [6]. The first one was by imposing an eigenvalue to be equal to +1 in the characteristic equation of the Jacobian of the discrete-time model. The second one is by expressing the steady-state feedback signal Fx(D) as a Fourier-series and then imposing the corresponding asymptotical critical condition. While the former approach is general and can be applied for PWM systems with A 1 A 2 , the last one has only been applied for systems with A 1 = A 2 . Its applicability to systems with A 1 A 2 has only been possible under some approximations. In both approaches, the external periodic function r(t) is assumed to be linear in such a way that its slope is constant during all the switching period. In this paper the steady-state approach is applied in time-domain for general systems with A 1 A 2 and B 1 B 2 and with a nonlinear external periodic function r(t) in such a way that the more general cases with a non constant slope within a switching period can also be treated without significant extra difficulty. Let us consider the following three different cases of compensating external signals.
where r 0 is the initial value of the external signal at t = 0 and σ e is its constant slope. This is the ideal case of most PWM systems. Since the slope is constant, the right-hand side of (18) is given by -Case of a quadratic modulating signal. For improving the performances of some switching PWM systems, a quadratic modulating signal can be used [21], [22]. Let σ 0 be the initial slope of the external signal at t = 0 and let σ T be its final slope at t = T . Therefore, this signal can be expressed as follows In this case, the slope is linearly dependent on the duty cycle D and the right-hand side of (18) is given by can be easily expressed as follows -In a practical implementation, the external modulating signal is implemented by an RC circuit making its shape more exponential than linear. In this case, the ramp signal r(t) can be expressed as follows where τ is a suitable time constant and σ 0 is the initial slope at the beginning of the switching cycle. The slope is exponentially depending on the duty cycle D and and the right-hand side of (18) is given by

Derivation of period-doubling (PD) bifurcation boundary in the time-domain
Consider a nonlinear dynamical system exhibiting a PD as shown in Fig. 3. After a PD bifurcation takes place, a 2T -periodic solution is born at the critical point while the T -periodic solution loses its stability but it still exists. During the switching cycle of duration T , a PWM system has two phases defined by the system matrices (A 1 , B 1 ) and (A 2 , B 2 ) respectively. During the switching cycle of duration 2T , a PWM system has four phases defined by the system matrices (A 1 , B 1 ), (A 2 , B 2 ), (A 1 , B 1 ) and (A 2 , B 2 ) respectively. Let us consider that the system is working at 2T -periodic orbit in steady-state. Therefore, during two consecutive switching periods in the interval (kT, (k+2)T ), let the crossing between the signals Fx(t) and r(t) occurs at t = (D − ε t + k)T and at t = (1 + D + ε t + k)T , k ∈ Z (see Fig. 4). The parameter ε t is a small quantity that vanishes at the boundary between T −periodic and 2T -periodic

PD bifurcation
Bifurcation parameter 2T -periodic T -periodic x Fig. 3. Sketch of a PD bifurcation in a nonlinear dynamical system and the corresponding waveforms before and after the bifurcation takes place by sweeping a parameter. behavior. At this point, the T -periodic solution and the 2Tperiodic solution are coincident (Fig. 3). By obtaining the expression of the 2T -periodic steady-state solutions at the switching instants, imposing the corresponding feedback constraints and equating them at the critical point (ε t → 0 ), a condition for predicting PD bifurcation is obtained in terms of the system matrices containing all the parameters. From the switching conditions at these two instants, the following equalities hold While the author in [6] and in some other papers, expresses the previous equations in the Fourier domain in the case (1−D)T 0 e A 2 t dtB 2 u and let Φ = Φ 1 Φ 2 . Exhibiting a 2Tperiodic regime, the sampled value of the steady-state state variables at the switching instants (D−ε t )T and (D+1+ε t )T can be obtained by using (3) and forcing 2T −periodicity. In doing so, it can be expressed as follows Subtracting (29) from (30), one obtains The boundary of PD bifurcation can be located by taking the limit in (40) when the parameter ε t → 0. Therefore at the onset of this instability the following equality holds While the right-hand side of (41) is generally easy to obtain, the left-hand side of the previous equation is mathematically more involved. Let us first focus on the righthand side of (41) and let us obtain the slope of the external signal for three different cases: -In the case of a linear ramp compensating signal, the right-hand side of (41) is given by -In the case of a quadratic modulating signal, the righthand side of (41) is given by can be easily expressed as follows (43) -In a practical implementation, the slope is exponentially depending on the duty cycle D and and the righthand side of (41) is given by As it was mentioned previously, the left-hand side of the previous equation is mathematically more involved. For simplicity let us consider that the external T −periodic function is linear during the switching period in such a way that its slope σ e is constant and that (41) can be written as follows where σ e,PD (D) the critical slope for PD bifurcation boundary is given by (46) By using (31)-(32), the limit expression in (46) becomes When the parameter ε t is very small (at the vicinity of the boundary of PD bifurcation), one has for i = 1, 2 in such a way that where the matrix M is given by the following expression Using the fact that (I±ε t Γ) −1 ≈ I∓ε t Γ for all small enough ε t and a suitable matrix Γ, the matrices Φ + (ε t ) and Φ − (ε t ) can still be simplified as follows It is worth noting here that the previous approximation is used only to calculate the limit expression as ε → 0 and that no approximation are used for calculating the rest of matrices. Therefore (47) becomes One can also demonstrate that in such a way that (54) becomes (57) Due to T periodicity when ε t → 0,one has that x((D + 1)T ) = x(DT ) and therefore (56) becomes where the subscript ss (for steady-state) in x has been omitted for simplicity. Also, one has that Φ 2 x(DT )+Ψ 2 = x(0). Therefore, rearranging terms and after some algebra, the critical slope σ e,PD (D) of the external periodic signal at the PD bifurcation boundary is given by where A = A 1 + A 2 and B = B 1 + B 2 . The previous expression can still be simplified after some arrangement using the expression of the matrix M and can be written as follows σ e,PD (D) = F[(I − Φ 2 ) −1 Φ 1 (I − Φ 2 Φ 1 )(f 1 (x(0)) + f 2 (x(0)))] (60) It is worth mentioning here that in [4], a slightly differently expressed condition has been obtained using a different approach based on solving the eigenvalue problem of the z-domain characteristic equation, for the same boundary condition which is reported and adapted here for comparison Again, although they are expressed differently, the critical ramp slope for PD bifurcation given in (60) and the one derived in [4] shown adapted in (61), are coincident.

Conclusions
In this work closed-form conditions for predicting the fastscale stability boundaries corresponding to both saddlenode (SN) and period-doubling (PD) bifurcation have been derived for a class of PWM switching systems with piecewise linear nonlinearity. Hence, the effect of the different parameters of the system upon the stability boundary can be unveiled. The general-purpose derived expressions can be applied to different examples of PWM systems such as switching power converters, switched capacitor chaos generators, temperature control systems and hydraulic valve drive control among others. The stability boundaries have been derived without the need of the Jacobian matrix and without expressing the system trajectories in the Fourier frequency-domain. The simple asymptotic time-domain approach used in this paper can be better understood by practitioners than those based in frequency-domain approach or on the eigenvalue problem of the Jacobian or the monodromy matrix.