# Stability by averaging via time-varying Lyapunov functions

Rami Katz, Frédéric Mazenc, Emilia Fridman

► To cite this version:

Rami Katz, Frédéric Mazenc, Emilia Fridman. Stability by averaging via time-varying Lyapunov functions. 22nd World Congress of the International Federation of Automatic Control (IFAC 2023), Jul 2023, Yokohama, Japan. 10.1016/j.ifacol.2023.10.1569. hal-04337744

HAL Id: hal-04337744

https://inria.hal.science/hal-04337744v1

Submitted on 13 Dec 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Stability by averaging via time-varying Lyapunov functions

Rami Katz**, Frédéric Mazenc* and Emilia Fridman**

* L2S-CNRS-CentraleSupélec, 3 rue Joliot Curie, 91192, Gif-sur-Yvette, France (email:frederic.mazenc@l2s.centralesupelec.fr)

** School of Electrical Engineering, Tel-Aviv University, Tel-Aviv (e-mail: rami@benis.co.il, emilia@eng.tau.ac.il)

**Abstract:** We study linear continuous-time systems with fast-varying almost periodic coefficients that are piecewise-continuous in time. Recently, a constructive time-delay approach to periodic averaging of systems with a single fast time-scale was introduced and employed to averaging of systems with small time-varying delays (of the order of the small parameter). In this paper we present a novel transformation of the fast varying coefficients. This transformation is suitable for averaging over multiple time-scales, and is applicable to averaging of systems with constant delays, where the value of delay is not small (i.e. essentially larger than the small parameter). We carry out stability analysis by employing time-varying Lyapunov functions (or functionals for the delayed case). The analysis leads to LMI conditions that are always feasible for small enough parameters. Numerical examples demonstrate the efficiency of the proposed approach and its conservatism.

**Keywords:** stability, averaging, time-varying.

## 1. INTRODUCTION

Systems with almost periodic coefficients and/or excitations are widely used in modern engineering applications (see e.g. Cheng et al. [2018], Sandberg and Möllerstedt [2001], Xie and Lam [2018] and the references therein). Such systems often include rapidly time-varying coefficients on multiple time-scales. Perturbation theory-based techniques have been useful in the analysis of systems with rapidly time-varying coefficients and have led to important qualitative results Bogoliubov and Mitropolskij [1961], Khalil [2001], Kokotovic and Khalil [1986], Teel et al. [2003], Moreau and Aeyels [2000].

The method of averaging has been especially efficient in the study of stability of systems with oscillatory control inputs Bullo [2002], Krstić and Wang [2000], Meerkov [1980]. The underlying idea behind asymptotic averaging is that stability of the original rapidly-varying system is guaranteed for small enough values of the parameter provided the averaged system is stable. However, by using such methods it is often difficult to find an explicit bound on the small parameter that preserves stability of the perturbed system. Direct Lyapunov techniques can be used to derive explicit bounds on the small parameter. For singularly perturbed systems, such bounds were derived in, e.g., Kokotovic and Khalil [1986] and Fridman [2002], via a direct Lyapunov approach.

Recently, a constructive time-delay approach to periodic averaging of a system with a single fast time-scale was introduced in Fridman and Zhang [2020]. It was suggested to use backward integration of the system, thereby pre-

senting it as a delayed neutral-type system, with the delay magnitude equal to the time-scale parameter. Stability of the delayed system then guarantees the stability of the original one. A direct Lyapunov method was used to obtain LMI-based conditions that allow to find an explicit upper bound on the small parameter which preserves the stability and ISS of the original system. This approach was also employed for averaging of systems with time-varying delays, where the delay bound which preserved the stability of the system was of the order of the time-scale parameter. These results were later extended to $L_2$-gain analysis for periodic averaging and to stochastic systems in Zhang and Fridman [2022]. However, the approach was confined to a single time-scale and the results for time-delay systems were limited to delays on the order of the small parameter.

In this paper we consider linear systems with rapidly-varying almost periodic coefficients that are piecewise-continuous in time. We begin by considering an ISS-like property of a system with two fast time-scales. Differently from Fridman and Zhang [2020], we employ a novel system transformation. This transformation leads to a new system whose ISS-like bounds guarantees the ISS-like bounds of the original one. For ISS-like analysis of the transformed system, we adopt a direct Lyapunov method leading to LMIs. In comparison to Fridman and Zhang [2020], our LMIs are significantly simpler. Feasibility of the LMIs provides quantitative estimates on the small parameters, which guarantee ISS-like estimates for the system. The LMIs are accompanied by theoretical guarantees on their feasibility for small values of the time-scale parameters.

We further extend our approach to rapidly time-varying systems with constant delay. Results on strict Lyapunov

<footnote>* Supported by Israel Science Foundation (grant no. 673/19) and by C. and H. Manderman Chair at Tel Aviv University</footnote>
functions for rapidly time-varying nonlinear systems were presented in Mazenc et al. [2006] and Mazenc and Malisoff [2017]. Differently from Fridman and Zhang [2020], the novel transformation suggested in this paper decouples the effects of the delay bound and time-scale parameter on the stability of the system, thereby allowing to obtain stability for non-small delay, relative to the time-scale parameter. We carry out stability analysis via a time-varying Lyapunov functional, leading to LMIs and provide theoretical guarantees on their feasibility.

**Notations:** Throughout the paper $\mathbb{R}^n$ denotes the n-dimensional Euclidean space with the vector norm $|\cdot|$, $\mathbb{R}^{n \times m}$ is the set of all $n \times m$ real matrices with the induced matrix norm $\|\cdot\|$. We also denote $\mathbb{Z}_+ = \{0, 1, 2, \dots\}$. The superscript $\top$ denotes matrix transposition, and the notation $P > 0$, for $P \in \mathbb{R}^{n \times n}$ means that $P$ is symmetric and positive definite. The symmetric elements of the symmetric matrix are denoted by $\ast$. For $0 < P \in \mathbb{R}^{n \times n}$ and $x \in \mathbb{R}^n$, we write $|x|^2_P = x^\top P x$. We denote by $W([-h, 0])$ the Banach space of a.e differentiable functions $\phi : [-h, 0] \to \mathbb{R}^n$ with square integrable derivative. The norm on $W([-h, 0])$ is given by the norm $\|\phi\|_W = \|\phi\|_W + \|\phi'\|_{L^2}$.

## 2. ISS OF RAPIDLY TIME-VARYING SYSTEMS

In this section we consider the fast-varying system

$$ \dot{x}(t) = A(t/\vartheta)x(t) + B(t/\epsilon)d(t), \quad t \ge 0. \quad (2.1) $$

Here $x(t) \in \mathbb{R}^n$ for $t \ge 0$, the small independent parameters $\vartheta, \epsilon > 0$ define two fast time-scales, $d \in C^1([0, \infty))$ is a disturbance and $A, B: \mathbb{R} \to \mathbb{R}^{n \times n}$ are uniformly bounded, piecewise continuous matrix functions (i.e. having at most a finite number of jump discontinuities on any interval in $\mathbb{R}$). The aim of this section is to derive explicit conditions, which guarantee ISS-like estimates for (2.1).

**Assumption 1:** $A$ and $B$ are almost periodic, meaning that there exist $0 < T$ and $A_{av}, B_{av}$, where $A_{av}$ is Hurwitz, such that for all $t \in \mathbb{R}$:

$$ T^{-1} \int_{t}^{t+T} L(\tau)d\tau = L_{av} + \Delta L(t), \quad L \in \{A, B\} \quad (2.2) $$

where $\Delta A, \Delta B: \mathbb{R} \to \mathbb{R}^{n \times n}$ are functions satisfying

$$ \sup_{\tau \in \mathbb{R}} \|\Delta B(\tau)\| \le \Delta_{B,M}, \quad \sup_{\tau \in \mathbb{R}} \|\Delta A(\tau)\| \le \Delta_{A,M} \quad (2.3) $$

for some $\Delta_{A,M}, \Delta_{B,M} > 0$.

**Assumption 2:** The time-scale parameters $\vartheta, \epsilon$ satisfy

$$ \max(\vartheta A_M, \epsilon B_M) \cdot T < 2 \quad (2.4) $$

where

$$ L_M := \sup_{\tau \in \mathbb{R}} \|L(\tau)\| < \infty, \quad L \in \{A, B\}. \quad (2.5) $$

Inspired by Mazenc et al. [2006], we introduce

$$ \begin{aligned} \omega_{\epsilon}(t) &= -\frac{1}{\epsilon T} \int_{t}^{t+\epsilon T} (t + \epsilon T - \tau) B(\tau/\epsilon) d\tau, \\ \varrho_{\vartheta}(t) &= -\frac{1}{\vartheta T} \int_{t}^{t+\vartheta T} (t + \vartheta T - \tau) A(\tau/\vartheta) d\tau, \quad t \ge 0. \end{aligned} \quad (2.6) $$

Note that piecewise continuity of $A$ and $B$ implies that $\omega_{\epsilon}(t)$ and $\varrho_{\vartheta}(t)$ are piecewise continuously differentiable. Indeed, both are continuously differentiable at any $t \ge 0$ for which $t/\epsilon + jT, j = 0, 1$ and $t/\vartheta + kT, k = 0, 1$ are

points of continuity of $B$ and $A$, respectively. Hence, any interval in $\mathbb{R}$ can be partitioned into finitely many subintervals with both $\omega_{\epsilon}(t)$ and $\varrho_{\vartheta}(t)$ piecewise continuously differentiable in the interior of each subinterval.

Using (2.5), (2.2) and (2.3), it can be easily verified that

$$ \sup_{\tau \in \mathbb{R}} \|\omega_{\epsilon}(\tau)\| \le 0.5\epsilon T B_M, \quad \sup_{\tau \in \mathbb{R}} \|\varrho_{\vartheta}(\tau)\| \le 0.5\vartheta T A_M \quad (2.7) $$

and for $t \ge 0$ (outside of a nowhere dense set):

$$ \begin{aligned} \dot{\omega}_{\epsilon}(t) &= B(t/\epsilon) - B_{av} - \Delta B(t/\epsilon), \\ \dot{\varrho}_{\vartheta}(t) &= A(t/\vartheta) - A_{av} - \Delta A(t/\vartheta), \end{aligned} \quad (2.8) $$

whereas by Assumption 2 and Neumann series, $I_n - \varrho_{\vartheta}(t)$ and $I_n - \omega_{\epsilon}(t)$ are invertible for all $t \in \mathbb{R}$ and

$$ \begin{aligned} \|(I_n - \varrho_{\vartheta}(t))^{-1}\| &= \left\| \sum_{k=0}^{\infty} \varrho_{\vartheta}(t)^k \right\|_{(2.7)} \le \frac{2}{2 - \vartheta T A_M}, \\ \|(I_n - \omega_{\epsilon}(t))^{-1}\| &= \left\| \sum_{k=0}^{\infty} \omega_{\epsilon}(t)^k \right\|_{(2.7)} \le \frac{2}{2 - \epsilon T B_M}. \end{aligned} \quad (2.9) $$

Using (2.6)-(2.8), the system (2.1) can be presented as

$$ \dot{x}(t) = [\dot{\varrho}_{\vartheta}(t) + A_{av} + \Delta A(t/\vartheta)] x(t) \\ + [\dot{\omega}_{\epsilon}(t) + B_{av} + \Delta B(t/\epsilon)] d(t), \quad t \ge 0. \quad (2.10) $$

To eliminate the terms $\dot{\varrho}_{\vartheta}(t)x(t)$ and $\dot{\omega}_{\epsilon}(t)d(t)$ in (2.10), we introduce the following transformation:

$$ z(t) = x(t) - \varrho_{\vartheta}(t)x(t) - \omega_{\epsilon}(t)d(t). \quad (2.11) $$

Note that $z(t)$ is continuous and piecewise continuously differentiable. In the proof of Theorem 2.1 we show that ISS-like estimates (with respect to $d$ and $\dot{d}$) on $z(t)$ imply ISS-like estimates for the system (2.1).

**Remark 2.1.** For the case of a single time-scale (i.e., $\epsilon = \vartheta$), differently from transformation (2.11), the time-delay transformation in Fridman and Zhang [2020] has a form

$$ z(t) = x(t) - G(t) \quad (2.12) $$

with

$$ G(t) = \frac{1}{\epsilon T} \int_{t-\epsilon T}^{t} (\tau - t + \epsilon T) [A(\tau)x(\epsilon\tau) + B(\tau)d(\epsilon\tau)] d\tau, $$

which leads to a neutral-type system. The transformation (2.12) allows for ISS of a larger class of essentially bounded disturbances than in the present paper.

Substituting (2.11) into (2.10), we obtain the following expression for $\dot{z}(t)$ (outside of a nowhere dense set):

$$ \dot{z}(t) = A_{av}z(t) + B_{av}d(t) + C_1(t)x(t) + C_2(t)d(t) - \omega_{\epsilon}(t)\dot{d}(t) \quad (2.13) $$

where

$$ \begin{aligned} C_1(t) &:= A_{av}\varrho_{\vartheta}(t) - \varrho_{\vartheta}(t)A(t/\vartheta) + \Delta A(t/\vartheta), \\ C_2(t) &:= A_{av}\omega_{\epsilon}(t) - \varrho_{\vartheta}(t)B(t/\epsilon) + \Delta B(t/\epsilon). \end{aligned} \quad (2.14) $$

Note that by (2.3) and (2.7), we have

$$ \begin{aligned} \sup_{\tau \in \mathbb{R}} \|C_1(\tau)\| &\le \frac{\vartheta T A_M [A_M + \|A_{av}\|]}{2} + \Delta_{A,M} := \bar{C}_1, \\ \sup_{\tau \in \mathbb{R}} \|C_2(\tau)\| &\le \frac{T B_M (\epsilon \|A_{av}\| + \vartheta A_M)}{2} + \Delta_{B,M} := \bar{C}_2. \end{aligned} \quad (2.15) $$

**Remark 2.2.** Differently from Fridman and Zhang [2020], where a fast-varying system over a single time-scale was considered and a neutral-type system was obtained by
integrating both sides of the ODE, here (2.11) allows us to treat the system (2.1), containing two time-scales. Furthermore, the approach presented in this manuscript allows us to avoid the introduction of a neutral-type system, which is more challenging for analysis.

For ISS-like analysis of (2.13), let $0 < \alpha$ and $0 < P \in \mathbb{R}^{n \times n}$. We introduce the Lyapunov function

$$V(t) = |z(t)|_P^2. \qquad (2.16)$$

Differentiating $V$ along the solution to (2.13), we have

$$\begin{aligned} \dot{V} + 2\alpha V &= z(t)^\top [PA_{av} + A_{av}^\top P + 2\alpha P] z(t) \\ &+ 2z(t)^\top PB_{av} d(t) + 2z(t)^\top PC_1(t)x(t) \\ &+ 2z(t)^\top PC_2(t)d(t) - 2z(t)^\top P\omega_\epsilon(t)\dot{d}(t). \end{aligned} \qquad (2.17)$$

Note that using Young's inequality, we have

$$\begin{aligned} |C_1(t)x(t)|^2 &\stackrel{(2.11)}{\le} \bar{C}_1^2 |(I - \varrho_\vartheta(t))^{-1} [z(t) + \omega_\epsilon(t)d(t)]|^2 \\ &\stackrel{(2.7)}{\le} \frac{2\bar{C}_1^2}{(2 - \vartheta T A_M)^2} [4|z(t)|^2 + (\epsilon T B_M)^2 |d(t)|^2], \\ |C_2(t)d(t)|^2 &\stackrel{(2.15)}{\le} \bar{C}_2^2 |d(t)|^2. \end{aligned} \qquad (2.18)$$

Let $0 < \lambda_i \in \mathbb{R}$, $i = 1, 2$. Taking into account (2.18), we apply the S-procedure by introducing

$$W_1(t) = \lambda_1 \left[ \frac{2\bar{C}_1^2}{(2 - \vartheta T A_M)^2} [4|z(t)|^2 + (\epsilon T B_M)^2 |d(t)|^2] - |C_1(t)x(t)|^2 \right] + \lambda_2 \left[ \bar{C}_2^2 |d(t)|^2 - |C_2(t)d(t)|^2 \right] \ge 0. \qquad (2.19)$$

Let $0 < \gamma_i \in \mathbb{R}$, $i = 1, 2$ and introduce for $t \ge 0$

$$\eta(t) = \operatorname{col} \{z(t), d(t), C_1(t)x(t), C_2(t)d(t), \omega_\epsilon(t)\dot{d}(t)\}. $$

From (2.17) and (2.19) we obtain

$$\begin{aligned} \dot{V} + 2\alpha V - \gamma_1^2 |d(t)|^2 - \gamma_2^2 |\omega_\epsilon(t)\dot{d}(t)|^2 + W_1 &\le \eta(t)^\top \Phi_{\epsilon, \vartheta} \eta(t) \le 0, \end{aligned} \qquad (2.20)$$

provided

$$\begin{aligned} \Phi_{\epsilon, \vartheta} &= \begin{bmatrix} \Phi_{\epsilon, \vartheta}^{(1)} & P & P & -P \\ * & 0 & 0 & 0 \\ * & -\lambda_1 I_n & 0 & 0 \\ * & * & -\lambda_2 I_n & 0 \\ * & * & * & -\gamma_2^2 I_n \end{bmatrix} < 0, \\ \Phi_{\epsilon, \vartheta}^{(1)} &= \begin{bmatrix} PA_{av} + A_{av}^\top P + 2\alpha P & PB_{av} \\ * & -\gamma_1^2 \end{bmatrix} + \Lambda_{\epsilon, \vartheta}, \qquad (2.21) \\ \Lambda_{\epsilon, \vartheta} &= \begin{bmatrix} \frac{8\lambda_1 \bar{C}_1^2}{(2 - \vartheta T A_M)^2} I_n & 0 \\ * & \lambda_2 \bar{C}_2^2 + \frac{2\lambda_1 (\epsilon T B_M \bar{C}_1)^2}{(2 - \vartheta T A_M)^2} \end{bmatrix}. \end{aligned}$$

Summarizing, we arrive at

**Theorem 2.1.** Consider the system (2.1) subject to Assumptions 1-2. Given the matrices $A_{av}, B_{av} \in \mathbb{R}^{n \times n}$ and constants $0 < \alpha, \epsilon, \vartheta, T, A_M, B_M, \Delta_{A,M}, \Delta_{B,M}$, let there exist $0 < P \in \mathbb{R}^{n \times n}$ and scalars $0 < \lambda_i, \gamma_i^2, i = 1, 2$, where $\lambda_i, i = 1, 2$ are given in (2.19), such that the LMI (2.21) holds. Then (2.1) satisfies the following ISS-like estimate:

$$\begin{aligned} |x(t)|^2 &\le \beta_1^2 e^{-2\alpha t} |x(0)|^2 + \beta_2^2 \max_{s \in [0,t]} |d(s)|^2 \\ &+ \beta_3^2 \max_{s \in [0,t]} |\dot{d}(s)|^2, \quad t \ge 0 \end{aligned} \qquad (2.22)$$

for some $\beta_i > 0$, $i = 1, 2, 3$. The LMI (2.21) is always feasible for small enough $0 < \alpha, \epsilon, \vartheta, \Delta_{A,M}, \Delta_{B,M}$ and large enough $0 < \gamma_i^2$, $i = 1, 2$. Feasibility of (2.21) for $\alpha, \epsilon, \vartheta, T, A_M, B_M, \Delta_{A,M}, \Delta_{B,M}, P, \lambda_i, \gamma_i^2, i = 1, 2$ implies its feasibility with any $\epsilon' < \epsilon, \vartheta' < \vartheta$ and the same $\alpha, T, A_M, B_M, \Delta_{A,M}, \Delta_{B,M}, P, \lambda_i, \gamma_i^2, i = 1, 2$, whence $\beta_i$ in (2.22) are independent of $\epsilon' \in (0, \epsilon)$ and $\vartheta' \in (0, \vartheta)$.

**Proof:** Fix $\tau > 0$. Feasibility of (2.21) implies that for all $t \in [0, \tau]$

$$\begin{aligned} \dot{V} + 2\alpha V - \gamma_1^2 |d(t)|^2 - \gamma_2^2 |\omega_\epsilon(t)\dot{d}(t)|^2 &\le 0 \\ \Rightarrow V(t) &\le e^{-2\alpha t} V(0) \\ &+ \int_0^t e^{-2\alpha(t-s)} \left( \gamma_1^2 |d(s)|^2 + \gamma_2^2 |\omega_\epsilon(s)\dot{d}(s)|^2 \right) ds. \end{aligned}$$

Since $\sigma_{min}(P) |z(t)|^2 \le V(t) \le \sigma_{max}(P) |z(t)|^2$ for all $t \ge 0$, we obtain for $t \ge 0$

$$|z(t)|^2 \stackrel{(2.7)}{\le} \frac{\sigma_{max}(P)}{\sigma_{min}(P)} e^{-2\alpha t} |z(0)|^2 + \frac{\gamma_1^2}{2\sigma_{min}(P)\alpha} \times \max_{s \in [0, \tau]} |d(s)|^2 + \frac{\gamma_2^2 (\epsilon T B_M)^2}{8\sigma_{min}(P)\alpha} \max_{s \in [0, \tau]} |\dot{d}(s)|^2, \qquad (2.23)$$

meaning that (2.13) satisfies an ISS-like estimate. Recalling (2.11), we have

$$\begin{aligned} |x(t)|^2 &\stackrel{(2.11)}{=} |(I - \varrho_\vartheta(t))^{-1} [z(t) + \omega_\epsilon(t)d(t)]|^2 \\ &\stackrel{(2.7)}{\le} \frac{8}{(2 - \vartheta T A_M)^2} |z(t)|^2 + \frac{2(\epsilon T B_M)^2}{(2 - \vartheta T A_M)^2} |d(t)|^2. \end{aligned} \qquad (2.24)$$

From (2.7) and (2.11), we find

$$|z(0)|^2 \le 2 \left(1 + \frac{\vartheta T A_M}{2}\right)^2 |x(0)|^2 + \frac{(\epsilon T B_M)^2}{2} |d(0)|^2. \qquad (2.25)$$

Finally, (2.23)-(2.25) imply (2.22), where $\beta_i, i = 1, 2, 3$ can be computed explicitly. We omit the formulas for $\{\beta_i\}_{i=1}^3$ due to space constraints.

Next, set $\epsilon = \vartheta = \Delta_{A,M} = \Delta_{B,M} = 0$ and consider $\Phi_1$ in (2.21) with $\Lambda_{\epsilon, \vartheta} = 0$ (see (2.19)). By Assumption 1, for $0 < \alpha$ small enough and $0 < \gamma_1^2$ large enough, we find that there exists $0 < P \in \mathbb{R}^{n \times n}$ such that $\Phi_1 < 0$. Fixing this $P$, applying Schur complement in $\Phi$ with respect to diag $\{-\lambda_1 I_n, -\lambda_2 I_n, -\gamma_2^2 I_n\}$ and taking $0 < \lambda_i, \gamma_i^2, i = 1, 2$ large enough we find that (2.21) holds. By continuity of eigenvalues, (2.21) remains true for the same values of $\alpha, P, \lambda_i, \gamma_i^2, i = 1, 2$ and small enough $0 < \epsilon, \vartheta, \Delta_{A,M}, \Delta_{B,M}$.

Assume (2.21) holds for $\alpha, \epsilon, \vartheta, T, A_M, B_M, \Delta_{A,M}, \Delta_{B,M}, P, \lambda_i, \gamma_i^2, i = 1, 2$. Replacing $\epsilon$ and $\vartheta$ by $\epsilon' < \epsilon$ and $\vartheta' < \vartheta$, respectively, and recalling (2.19), we have $\Lambda_{\epsilon', \vartheta'} \le \Lambda_{\epsilon, \vartheta}$. Hence, we obtain $\Phi_{\epsilon', \vartheta'}^{(1)} \le \Phi_{\epsilon, \vartheta}^{(1)}$ in (2.21). Taking the latter into account and applying Schur complement to $\Phi_{\epsilon', \vartheta'}$ in (2.21) we find that $\Phi_{\epsilon, \vartheta} < 0$ implies $\Phi_{\epsilon', \vartheta'} < 0$. $\square$

### 3. RAPIDLY TIME-VARYING SYSTEMS WITH CONSTANT DELAY

We consider the system
$$
\begin{equation}
\begin{aligned}
\dot{x}(t) &= Ax(t) + B(t/\epsilon)x(t-h), & t \ge 0, \\
x(\theta) &= \phi(\theta), & \theta \in [-h, 0]
\end{aligned}
\tag{3.1}
\end{equation}
$$

where $x(t) \in \mathbb{R}^n$ for $t \ge 0$, $A \in \mathbb{R}^{n \times n}$, $h, \epsilon > 0$ and $\phi \in W([-h, 0])$. Here, the delayed term $x(t-h)$ is multiplied by a fast-varying coefficient $B(t/\epsilon)$, where $B: \mathbb{R} \to \mathbb{R}^{n \times n}$ is a piecewise continuous function subject to Assumption 1 and inequality (2.5). We further make the following assumption

**Assumption 3:** The matrix $A_0 = A + B_{av}$ is Hurwitz.

Remark 3.1. We assume that A in (3.1) is constant for simplicity of presentation only. Our approach can be easily extended to the case of a system of the form

$$
\dot{x}(t) = A(t/\vartheta)x(t) + B(t/\epsilon)x(t-h), \quad t \ge 0
$$

where $\epsilon, \vartheta > 0$ and $A : \mathbb{R} \to \mathbb{R}^{n \times n}$ is a piecewise
continuous function satisfying Assumption 1 and $A_{av}+B_{av}$
is Hurwitz. In this case, the transformation (3.3) below will
be replaced by the following transformation:

$$
z(t) = (I - \varrho_{\vartheta}(t)) x(t) - \omega_{\epsilon}(t)x(t - h), \quad t \ge 0
$$

with $\varrho_{\vartheta}(t)$ and $\omega_{\epsilon}(t)$ given in (2.6). In this case $A_{av}$ need not by Hurwitz.

Using $\omega_{\epsilon}(t)$ in (2.6), $\dot{\omega}_{\epsilon}(t)$ in (2.8) and (2.2), (3.1) can be
presented as

$$
\dot{x}(t) = Ax(t) + [\dot{\omega}_{\epsilon}(t) + B_{av} + \Delta B (t/\epsilon)] x(t-h) \quad (3.2)
$$

Introduce

$$
z(t) = x(t) - \omega_{\epsilon}(t)x(t-h), \quad t \ge 0 \quad (3.3)
$$

which is again continuous and piecewise continuously
differentiable. Differentiating $z(t)$ for $t \ge h$ (outside of
a nowhere dense set), we obtain the following:

$$
\begin{equation}
\begin{aligned}
\dot{z}(t) ={}& A_0 z(t) - B_{av} \zeta_h(t) + \Delta B (t/\epsilon) x(t-h) \\
& + H_\epsilon(t) x(t-h) + G_\epsilon(t) x(t-2h), \quad t \ge h
\end{aligned}
\tag{3.4}
\end{equation}
$$

where $A_0$ is given in Assumption 3,

$$
\zeta_h(t) = z(t) - z(t-h), \quad t \ge h, \qquad (3.5)
$$

is state error induced by the delay and

$$
\begin{align}
H_{\epsilon}(t) &= A\omega_{\epsilon}(t) - \omega_{\epsilon}(t)A, \nonumber \\
G_{\epsilon}(t) &= B_{av}\omega_{\epsilon}(t-h) - \omega_{\epsilon}(t)B((t-h)/\epsilon). \tag{3.6}
\end{align}
$$

By using (2.6), (2.7) and (2.8), we see that

$$
\begin{align}
\sup_{\tau \in \mathbb{R}} \| H_{\epsilon}(\tau) \| &\le \epsilon T B_{M} \| A \| := \bar{H}_{\epsilon}, \nonumber \\
\sup_{\tau \in \mathbb{R}} \| G_{\epsilon}(\tau) \| &\le 0.5 \epsilon T B_{M} (\| B_{av} \| + B_{M}) :=: \bar{G}_{\epsilon}. \tag{3.7}
\end{align}
$$

For exponential stability of $z(t)$ subject to (3.4), let $0 < P, S, R \in \mathbb{R}^{n \times n}$ and $0 < \alpha, q_i, i = 1, \dots, 4$. We consider the Lyapunov functional

$$
V(t) = |z(t)|_{P}^{2} + V_{S}(t) + V_{R}(t) + \sum_{i=1}^{4} V_{q_{i}}(t),
$$

where

$$
\begin{align*}
V_S(t) &= \int_{t-h}^{t} e^{-2\alpha(t-s)} |z(s)|_S^2 ds, \\
V_R(t) &= h \int_{-h}^{0} \int_{t+\theta}^{t} e^{-2\alpha(t-s)} |\dot{z}(s)|_R^2 ds d\theta,
\end{align*}
$$

are introduced to compensate $\zeta_h(t)$ and

$$
\begin{align}
V_{q_1}(t) &= q_1 \int_{t-2h}^{t} e^{-2\alpha(t-s)} |G_{\epsilon}(s+2h)x(s)|^2 ds, \nonumber \\
V_{q_2}(t) &= q_2 \int_{t-h}^{t} e^{-2\alpha(t-s)} |\omega_{\epsilon}(s+h)x(s)|^2 ds, \nonumber \\
V_{q_3}(t) &= q_3 \int_{t-h}^{t} e^{-2\alpha(t-s)} |H_{\epsilon}(s+h)x(s)|^2 ds, \nonumber \\
V_{q_4}(t) &= q_4 \int_{t-h}^{t} e^{-2\alpha(t-s)} |\Delta B((s+h)/\epsilon)x(s)|^2 ds \tag{3.8}
\end{align}
$$

will compensate $\Delta B(\frac{t}{\epsilon}) x(t-h)$, $H_{\epsilon}(t)x(t-h)$ and $G_{\epsilon}(t)x(t-2h)$.

Differentiating $V_S$ and $V_R$ along the solution to (3.4) gives

$$
\begin{equation}
\begin{split}
\dot{V}_S + 2\alpha V_S ={}& (1 - e^{-2\alpha h}) |z(t)|_S^2 - e^{-2\alpha h} |\zeta_h(t)|_S^2 \\
& + 2e^{-2\alpha h} z(t)^\top S \zeta_h(t),
\end{split}
\tag{3.9}
\end{equation}
$$

$$
\dot{V}_R + 2\alpha V_R \le h^2 |\dot{z}(t)|_R^2 - e^{-2\alpha h} |\zeta_h(t)|_R^2
$$

whereas

$$
\begin{equation}
\begin{split}
\frac{d}{dt} |z(t)|_P^2 + 2\alpha |z(t)|_P^2 ={}& z(t)^\top [P A_0 + A_0^\top P + 2\alpha P] z(t) \\
& - 2z(t)^\top P B_{av} \zeta_h(t) + 2z(t)^\top P \Delta B (t/\epsilon) x(t-h) \\
& + 2z(t)^\top P H_{\epsilon}(t) x(t-h) + 2z(t)^\top P G_{\epsilon}(t) x(t-2h).
\end{split}
\tag{3.10}
\end{equation}
$$

Differentiation of $V_{q_1}(t)$ gives

$$
\begin{equation}
\begin{split}
\dot{V}_{q_1} + 2\alpha V_{q_1} ={}& q_1 |G_{\epsilon}(t+2h)x(t)|^2 \\
& - q_1 e^{-4\alpha h} |G_{\epsilon}(t)x(t-2h)|^2 \stackrel{(3.3)}{=} q_1 \bar{G}_{\epsilon}^2 |z(t)|^2 \\
& + q_1 \bar{G}_{\epsilon}^2 |\omega_{\epsilon}(t)x(t-h)|^2 + 2q_1 \bar{G}_{\epsilon}^2 \times \\
& z(t)^{\top} \omega_{\epsilon}(t)x(t-h) - q_1 e^{-4\alpha h} |G_{\epsilon}(t)x(t-2h)|^2
\end{split}
\tag{3.11}
\end{equation}
$$

whereas similar arguments for $V_{q_i}$, $i = 2, 3, 4$ yield

$$
\begin{equation}
\begin{aligned}
\dot{V}_{q_2} + 2\alpha V_{q_2} \le & q_2 \frac{(\epsilon T B_M)^2}{4} |z(t)|^2 + q_2 \frac{(\epsilon T B_M)^2}{2} z(t)^\top \times \\
& \omega_\epsilon(t) x(t-h) - q_2 \left[ e^{-2\alpha h} - \frac{(\epsilon T B_M)^2}{4} \right] |\omega_\epsilon(t) x(t-h)|^2, \\
\dot{V}_{q_3} + 2\alpha V_{q_3} \le & q_3 \bar{H}_\epsilon^2 |z(t)|^2 + q_3 \bar{H}_\epsilon^2 |\omega_\epsilon(t)x(t-h)|^2 \\
& + 2q_3 \bar{H}_\epsilon^2 z(t)^\top \omega_\epsilon(t)x(t-h) - q_3 e^{-2\alpha h} |H_\epsilon(t)x(t-h)|^2, \\
\dot{V}_{q_4} + 2\alpha V_{q_4} \le & q_4 \Delta_{B,M}^2 |z(t)|^2 + q_4 \Delta_{B,M}^2 |\omega_\epsilon(t)x(t-h)|^2 \\
& + 2q_4 \Delta_{B,M}^2 z(t)^\top \omega_\epsilon(t)x(t-h) - q_4 e^{-2\alpha h} \times \\
& |\Delta B (t/\epsilon) x(t-h)|^2.
\end{aligned}
\tag{3.12}
\end{equation}
$$

Finally,

$$
\begin{align}
h^2 |\dot{z}(t)|_{R_z}^2 &= h^2 \eta^\top(t) \Lambda^\top R \Lambda \eta(t), \nonumber \\
\Lambda &= [A_0, -B_{av}, I_n, I_n, I_n, 0], \nonumber \\
\eta(t) &= \operatorname{col}\{z(t), \zeta_h(t), \Delta B(t/\epsilon)x(t-h), H_\epsilon(t)x(t-h) \nonumber \\
& \phantom{=} G_\epsilon(t)x(t-2h), \omega_\epsilon(t)x(t-h)\}. \tag{3.13}
\end{align}
$$

From (3.10)-(3.13) we have

$$
\dot{V} + 2\alpha V \le \eta(t)^{\top} \Psi_{\epsilon,h} \eta(t) \le 0 \quad (3.14)
$$

provided
$$
\Psi_{\epsilon,h} = \begin{bmatrix} \Psi_{\epsilon}^{(1)} & P & P & P & \sigma_{\epsilon} I_{n} \\ 0 & 0 & 0 & 0 & 0 \\ * & -q_{4}e^{-2\alpha h} I_{n} & 0 & 0 & 0 \\ * & * & -q_{3}e^{-2\alpha h} I_{n} & 0 & 0 \\ * & * & * & -q_{1}e^{-4\alpha h} I_{n} & 0 \\ * & * & * & * & -\beta_{\epsilon} I_{n} \\ \end{bmatrix} + h^{2}\Lambda^{\top}R\Lambda < 0,
$$

$$
\Psi_{\epsilon}^{(1)} = \begin{bmatrix} \psi & -PB_{av} + e^{-2\alpha h}S \\ * & -e^{-2\alpha h}(S+R) \end{bmatrix} + \sigma_{\epsilon} \begin{bmatrix} I_{n} & 0 \\ 0 & 0 \end{bmatrix},
$$

$$
\psi = PA_0 + A_0^\top P + 2\alpha P + (1 - e^{-2\alpha h}) S,
$$

$$
\sigma_{\epsilon} = \frac{q_2(\epsilon TB_M)^2}{4} + q_1 \bar{G}_{\epsilon}^2 + q_3 \bar{H}_{\epsilon}^2 + q_4 \Delta_{B,M}^2,
$$

$$
\beta_{\epsilon} = q_2 \frac{4e^{-2\alpha h} - (\epsilon TB_M)^2}{4} - q_1 \bar{G}_{\epsilon}^2 - q_3 \bar{H}_{\epsilon}^2 - q_4 \Delta_{B,M}^2. \quad (3.15)
$$

We are now ready to state the main theorem of this section.
Theorem 3.1. Consider the system (3.1) subject to Assumptions 1 and 3. Given the matrices $A_0, A, B_{av} \in \mathbb{R}^{n \times n}$, constants $0 < \alpha, \epsilon, T, B_M, h$ and tuning parameter $0 < \Delta_{B,M}$, let there exist matrices $0 < P, S, R \in \mathbb{R}^{n \times n}$ and scalars $0 < q_i, i = 1, \dots, 4$ such that $\kappa_{\epsilon} := \frac{\epsilon TB_M}{2} < e^{-\alpha h}$ and the LMI (3.15) hold. Then the system (3.1) is exponentially stable with decay rate $\alpha > 0$. The LMI (3.15) is always feasible for small enough $0 < \alpha, \epsilon, h, \Delta_{B,M}$. Furthermore, let (3.15) hold for some $\epsilon, h, \alpha, T, B_M, \Delta_{B,M}, P, S, R, q_i, i = 1, \dots, 4$. Then, (3.15) is feasible for any $\bar{\epsilon} < \epsilon$ and the same $h, \alpha, T, B_M, \Delta_{B,M}, P, S, R, q_i, i = 1, \dots, 4$.

**Proof:** Due to space constraints, we only provide a sketch.
First, Feasibility of the LMI (3.15) implies that there exists
a constant $0 < M$ such that

$$
|z(t)| \leq Me^{-\alpha(t-h)}, \quad t \geq h. \tag{3.16}
$$

Next, for any $k \in \mathbb{N}$, we denote $X_k = \sup_{\tau \in [kh,(k+1)h]} |x(\tau)|$.
From (3.3) and (3.16), we find that $X_{k+1} \le Me^{-\alpha kh} + \kappa_{\epsilon}X_k$, $k \in \mathbb{N}$. Comparing with the solution of the linear
difference equation

$$
Y_{k+1} = Me^{-\alpha kh} + \kappa_{\epsilon}Y_k, \quad k \in \mathbb{N}. \quad (3.17)
$$

with $X_1 = Y_1$ it follows that

$$
|x(t)| \leq X_k \leq (\kappa_{\epsilon}^{-1} (X_1 - \mu_{\epsilon,h}) + \mu_{\epsilon,h} e^{\alpha h}) e^{-\alpha(t-h)}.
$$

Consider (3.15) with $\epsilon = \Delta_{B,M} = 0$, which implies $\sigma_{\epsilon} = 0$ and $\beta_{\epsilon} = q_{2}e^{-2\alpha h}$. By Assumption 3, there exist $0 < P, S, R \in \mathbb{R}^{n \times n}$ and small $0 < \alpha \in \mathbb{R}$ such that

$$
\Psi_1 + h^2 \begin{bmatrix} A_0^\top \\ -B_{av}^\top \end{bmatrix} R [A_0 - B_{av}] < 0, \quad (3.18)
$$

provided $h > 0$ is small enough. Fixing such $0 < \alpha, h$,
taking $q_i = q$, $i = 1, \dots, 4$ for $q > 0$ large enough
and applying Schur complement in $\Psi_{\epsilon,h}$, we find that
(3.15) holds. By continuity of eigenvalues, we have that
$\kappa_{\epsilon} < e^{-\alpha h}$ and (3.15) hold for the same $\alpha, h$ and small
enough $0 < \epsilon, \Delta_{B,M}$.

If (3.15) hold with $\epsilon, h, \alpha, T, B_M, \Delta_{B,M}, P, S, R, q_i, i = 1, \dots, 4$ then Schur complement of (3.15) and monotonicity of $\kappa_{\epsilon} < e^{-\alpha h}$ yield feasibility of (3.15) with $\bar{\epsilon} < \epsilon$ and the same $h, T, B_M, \Delta_{B,M}, P, S, R, q_i, i = 1, \dots, 4$. $\square$

**Corollary 1.** Fix $\alpha, T, B_M, \Delta_{B,M}$ and let $\epsilon_* > 0$ be such that there exists $h > 0$ for which the LMI (3.15) is feasible

with $(\epsilon_*, h)$. Define $h_* : (0, \epsilon_*) \to (0, \infty)$ by

$$
h_{*}(\epsilon) = \sup \{h > 0 \mid (3.15) \text{ is feasible for } (\epsilon, h)\}. \quad (3.19)
$$

Then $h_*(\epsilon)$ is positive and non-increasing on $(0, \epsilon_*)$.

Remark 3.2. In Fridman and Zhang [2020], the delay bound which preserved the stability of the system was of the order of the time-scale $\epsilon$. Here (3.3) decouples the effects of the delay bound and $\epsilon$ on the stability of the system, thereby allowing to obtain stability of (3.1) for non-small delay $h$, relative to $\epsilon$.

4. NUMERICAL EXAMPLE

We show that our method is complementary to that of
Fridman and Zhang [2020] and Zhang and Fridman [2022].

4.1 Rapidly varying persistently excited system

We consider stability of the following persistently excited
system [Fridman and Zhang 2020, Example 3.1]

$$
\dot{x}(t) = -p(t/\vartheta)p(t/\vartheta)^{\top}x(t), \quad t \ge 0,
$$

where $x(t) \in \mathbb{R}^n$ and there exists constants $0 < \rho < M$
such that

$$
\rho I_n \le \int_0^1 p(\tau - \theta) p(\tau - \theta)^\top \le M^2 I_n, \quad \tau \ge 1.
$$

Recalling (2.1), here we have $A(t/\vartheta) = -p(t/\vartheta)p(t/\vartheta)^{\top}$,
whereas $B \equiv 0$.

We consider the case $n = 1$, $A_{av} = 0.5(M^2 + \rho)$, $\Delta A = 0.5(M^2 - \rho)$, where $M = 1$ and $\rho = 0.55$, as in [Fridman and Zhang 2020, Example 3.1]. Choosing further $\alpha = 0.5$, we find that feasibility of the LMI of Theorem 2.1 is preserved for $\vartheta_{max} = 0.1544$. In comparison to [Fridman and Zhang 2020, Example 3.1], where $\vartheta_{max} = 0.0645$, our approach performs better in this example.

4.2 Stabilization by fast switching

Consider stabilization by fast switching of a linear system
[Fridman and Zhang 2020, Example 4.2]. Let $\epsilon > 0$ and

$$
A_1 = \begin{bmatrix} 0.1 & 0.3 \\ 0.6 & -0.2 \end{bmatrix}, \quad A_2 = \begin{bmatrix} -0.13 & -0.16 \\ -0.33 & 0.03 \end{bmatrix}
$$

and set

$$
A(\tau) = \begin{cases} A_1, & \tau \in [k, k+0.4), \ k \in \mathbb{Z}_+ \\ A_2, & \tau \in [k+0.4, k+1), \ k \in \mathbb{Z}_+ \end{cases} . \quad (4.1)
$$

Given $\tau \in [k, k+1)$, (4.1) can be written as

$$
A(\tau) = \chi_{[k, k+0.4)}(\tau) A_1 + [1 - \chi_{[k+0.4, k+1)}(\tau)] A_2,
$$

where $\chi_{[a,b]}$ is the indicator function of the interval $[a,b)$.
Note that $A(\tau)$ is 1-periodic.

Consider (2.1) with $\vartheta = \epsilon$, $A(t/\vartheta)$ from (4.1), $B(t) = [0, 1]^\top$. Recalling (2.2) and (2.5), we find that $T = 1$ and

$$
A_{av} = 0.4A_1 + 0.6A_2 = \begin{bmatrix} -0.038 & 0.024 \\ 0.042 & -0.062 \end{bmatrix}, \quad B_{av} = B
$$

where $\sigma_{max}(A_{av}) \approx -0.08394$. Furthermore, $A_M = \max(\|A_1\|, \|A_2\|) \approx 0.63245$ and $B_M = 1$. We fix $\alpha = 0.005$ and verify the LMI of Theorem 2.1. The LMI was found feasible for $\vartheta_{max} = 0.0376$. Next, we consider $\epsilon = 0.025$ constant disturbance $d(t) \equiv d$. An upper bound
<figure><img src="image_2.png" /><figcaption>Fig. 1. Theorem 3.1: feasibility preserving $(\epsilon, h)$. Top: $X = h, Y = \max \epsilon$. Bottom: $X = \epsilon, Y = \max h$</figcaption></figure>

on the ISS gain, $\beta_2$ in (2.22) (computed from (2.23)-(2.25)), is given by $1.079 \times 10^2$. In comparison to [Zhang and Fridman 2022, Example 4.2], where the corresponding value of $\vartheta_{max} = 0.1306$, the results obtained using the approach developed in this paper are conservative. However, differently from Fridman and Zhang [2020], our method can easily cope with two fast time-scales.

Consider now the system (2.1) with a constant disturbance $d(t) \equiv d$, $\epsilon \neq \vartheta$, $A(\tau)$ from (4.1) and $B(\theta) = [0 \ 1 + 0.1 \sin(2\pi\theta)]^\top$. Note that $B(\theta)$ is 1-periodic with $B_{av} = [0 \ 1]^\top$ and $B_M = 1.1$. We fix $\vartheta = 0.025$ and $\epsilon \in \frac{1}{\pi} \cdot [0.0002, 0.5]$ with grid of step side 0.0001. Note that $A(t/\vartheta)$ and $B(t/\epsilon)$ do not have a common period. Hence, the approach in Fridman and Zhang [2020] cannot be used in this system. Choosing $\alpha = 0.005$ and restricting the bound on the ISS gain, $\beta_2$ in (2.22) (computed from (2.23)-(2.25)), to be at most $1.079 \times 10^2$, we verify the LMI of Theorem 2.1 to obtain the maximal value of $\epsilon$ which preserves feasibility of the LMI. The LMI was found feasible for $\epsilon_{max} = 0.178$. Comparing to the previous simulation with $\vartheta_{max} = \epsilon_{max} = 0.0376$, we see that for two independent time-scale parameters, the magnitude of $\vartheta$ is a larger constraint for feasibility of the LMI of Theorem 2.1 than the magnitude of $\epsilon$ in this example.

Next, we consider the delayed system (3.1), where we choose $A = 0$ and $B(t) = A(t)$ from (4.1). We fix the decay rate $\alpha = 0.005$. Then, the LMIs of Theorem 3.1 were checked for feasibility twice. First, for $\epsilon \in [0, 0.0475]$ the LMIs were verified to find the maximal value of $h$, which preserves feasibility. Second, for $h \in [0, 0.356]$ the LMIs were verified to find the maximal value of $\epsilon$, which preserves feasibility. The results are given in Figure 1. Differently from Fridman and Zhang [2020], where $\epsilon$ and the delay bound were of the same order, it can be seen that as $\epsilon \to 0^+$, the value of $h_*(\epsilon)$ (see (3.19)) *grows*.

## 5. CONCLUSION

This paper presented a novel quantitative approach to stability of linear continuous-time systems. Differently from the time-delay approach to averaging, our method is efficient for systems with fast-varying almost periodic coefficients with multiple independent time-scales, as well as to such systems in the presence of constant delay, which does not vanish when the small parameter approaches zero. Moreover, our approach leads to essentially simpler single reduced-order LMI with fewer decision parameters. Future work may include improvement of the method and its extension to time-varying delays.

## REFERENCES

* Bogoliubov, N.N. and Mitropolskij, I.A. (1961). *Asymp-totic methods in the theory of non-linear oscillations*, volume 10. CRC Press.

* Bullo, F. (2002). Averaging and vibrational control of mechanical systems. *SIAM Journal on Control and Optimization*, 41(2), 542–562.

* Cheng, X., Tan, Y., and Mareels, I. (2018). On robustness analysis of linear vibrational control systems. *Automatica*, 87, 202–209.

* Fridman, E. (2002). Effects of small delays on stability of singularly perturbed systems. *Automatica*, 38(5), 897–902.

* Fridman, E. and Zhang, J. (2020). Averaging of linear systems with almost periodic coefficients: A time-delay approach. *Automatica*, 122, 109287.

* Khalil, H.K. (2001). Nonlinear systems third edition.

* Kokotovic, P.V. and Khalil, H.K. (1986). *Singular pertur-bations in systems and control*. IEEE press.

* Krstić, M. and Wang, H.H. (2000). Stability of extremum seeking feedback for general nonlinear dynamic systems. *Automatica*, 36(4), 595–601.

* Mazenc, F. and Malisoff, M. (2017). Extensions of Razumikhin's theorem and Lyapunov-Krasovskii functional constructions for time-varying systems with delay. *Automatica*, 78, 1–13.

* Mazenc, F., Malisoff, M., and De Queiroz, M.S. (2006). Further results on strict Lyapunov functions for rapidly time-varying nonlinear systems. *Automatica*, 42(10), 1663–1671.

* Meerkov, S. (1980). Principle of vibrational control: theory and applications. *IEEE Transactions on Automatic Control*, 25(4), 755–762.

* Moreau, L. and Aeyels, D. (2000). Practical stability and stabilization. *IEEE Transactions on Automatic Control*, 45(8), 1554–1558.

* Sandberg, H. and Möllerstedt, E. (2001). Periodic modelling of power systems. *IFAC Proceedings Volumes*, 34(12), 89–94.

* Teel, A.R., Moreau, L., and Nesic, D. (2003). A unified framework for input-to-state stability in systems with two time scales. *IEEE Transactions on Automatic Control*, 48(9), 1526–1544.

* Xie, X. and Lam, J. (2018). Guaranteed cost control of periodic piecewise linear time-delay systems. *Automatica*, 94, 274–282.

* Zhang, J. and Fridman, E. (2022). $L_2$-gain analysis via time-delay approach to periodic averaging with stochastic extension. *Automatica*, 137, 110126.