Transfer Operator Based Approach for Optimal Stabilization of Stochastic System

TL;DR

This paper introduces a transfer operator-based method for the optimal stabilization of stochastic nonlinear systems, extending deterministic approaches with minimal modifications and demonstrating effectiveness through simulations.

Contribution

It develops a transfer Perron Frobenius operator approach for stochastic systems, bridging deterministic and stochastic stabilization techniques with a new computational framework.

Findings

Framework successfully stabilizes stochastic nonlinear systems

Finite-dimensional approximations enable practical computations

Simulation results confirm the effectiveness of the approach

Abstract

In this paper we develop linear transfer Perron Frobenius operator-based approach for optimal stabilization of stochastic nonlinear system. One of the main highlight of the proposed transfer operator based approach is that both the theory and computational framework developed for the optimal stabilization of deterministic dynamical system in [1] carries over to the stochastic case with little change. The optimal stabilization problem is formulated as an infinite dimensional linear program. Set oriented numerical methods are proposed for the finite dimensional approximation of the transfer operator and the controller. Simulation results are presented to verify the developed framework.

Equations

$$x_{n+1}=T(x_n,{\xi }_n),$$

$$\mathrm{Prob}({\xi }_n\in B)=v(B),\forall n,B\subset W,$$

$$[{\mathbb{P}}_T\mu ](A)={\int }_X\left\{{\int }_W{\chi }_A(T(x,y))dv(y)\right\}d\mu (x)$$

$$\mathrm{Prob}\{T^n(x_0,{\xi }_0^n)\in \mathcal{O}\}=1,\forall n\ge 0,$$

$$\mathrm{Prob}\{\underset{n\to \infty }{\lim }{T}^n(x_0,{\xi }_0^n)=0\}=1.$$

$$m\{x\in X:Prob\{T^n(x,{\xi }_0^n)\in B\}\ge {\alpha }^n\}\le K{\beta }^n,$$

$$[{\mathbb{P}}_T\mu ](B)$$

$$[{\mathbb{P}}_T^1\mu ](B)={\int }_{E^c}{\int }_W{\chi }_B(T(x,y))dv(y)d\mu (x),$$

$$[{\mathbb{P}}_T^0\mu ](B)={\int }_B{\chi }_B(T(x,y))dv(y)d\mu (x),$$

$$\mathbb{P}_{T}=\left[\begin{array}[]{cc}\mathbb{P}_{T}^{0}&0\\ \times&\mathbb{P}_{T}^{1}\end{array}\right].$$

$$[{\mathbb{P}}_T^1\bar{\mu }](B)<\gamma \bar{\mu }(B)$$

$$\gamma [{\mathbb{P}}_T^1\bar{\mu }](B)-\bar{\mu }(B)=-m(B)$$

$$x_{n+1}=T(x_n,u_n,{\xi }_n)=:T_{{\xi }_n}(x_n,u_n)$$

$$\theta (D):=[{\mathbb{P}}_C\mu ](D)=\mu (C^{-1}(D))$$

$$[{\mathbb{P}}_{C^{-1}}\theta ](B):=\mu (B)=\theta (C(B)),$$

$$x_{n+1}=T(x_n,K(x_n),{\xi }_n)=T(C(x_n),{\xi }_n)=:T_{{\xi }_n}\circ C(x_n).$$

$$[{\mathbb{P}}_T\theta ](B)={\int }_Y\left\{{\int }_W{\chi }_B(T(y,\xi ))dv(\xi )\right\}d\theta (y).$$

$$[{\mathbb{P}}_{T_{\xi }\circ C}\mu ](B)={\int }_Y\left\{{\int }_W{\chi }_B(T(y,\xi ))dv(\xi )\right\}d[{\mathbb{P}}_C\mu ](y)$$

$$[{\mathbb{P}}_T{\mathbb{P}}_C\mu ](B)={\int }_X{\int }_{C(x)}\left\{{\int }_W{\chi }_B(T(y,\xi ))dv(\xi )\right\}d{\theta }_x(y)d\mu (x).$$

$$C(x)=\left\{\begin{array}[]{ccl}C_{0}(x)=(x,K_{0}(x))&{\rm for}&x\in U(\epsilon)\\ C_{1}(x)=(x,K_{1}(x))&{\rm for}&x\in X_{1}.\end{array}\right.$$

$${\mathcal{C}}_C(B)={\int }_B\sum _{n=0}^{\infty }{\gamma }^n{E}_{{\xi }_0^n}[G(C(x_n),{\xi }_n)]dm(x),$$

$${\mathcal{C}}_C(B)={\int }_B\sum _{n=0}^{\infty }{\gamma }^n{E}_{{\xi }_0^n}[G(C(x_n),{\xi }_n)]dm(x)$$

$$={\int }_{\mathcal{Z}}G(y,\xi )d[{\mathbb{P}}_C^1{\bar{\mu }}_B](y)dv(\xi )={⟨G,d{\mathbb{P}}_C^1{\bar{\mu }}_{B}dv⟩}_{\mathcal{Z}},$$

$$\gamma [{\mathbb{P}}_T^1\cdot {\mathbb{P}}_C^1{\bar{\mu }}_B](D)-{\bar{\mu }}_B(D)=-m_B(D),$$

$${\mathcal{C}}^{\ast }(B)=\underset{C}{\min }{\mathcal{C}}_C(B).$$

$$\underset{\theta \ge 0}{\min }$$

$$\gamma [{\mathbb{P}}_T^1\theta ](D)-[{\mathbb{P}}_{P_1}^1\theta ](D)=-m_B(D),$$

$$Prob\{{\xi }_n={\xi }^{\ell }\}=v^{\ell },\forall n,\ell =1,\dots ,L.$$

$$(P_{T_{a,\ell }}{)}_{(ij)}:=\frac{m(T_{u^a,{\xi }^{\ell }}^{-1}(D_j)\cap D_i)}{m(D_i)},$$

$$\underset{{\theta }^a,\mu \ge 0}{\min }$$

$$\gamma \sum _{a=1}^M\left[\sum _{\ell =1}^L{v}^{\ell }(P_{T_{a,\ell }}{)}^{{}^{\prime }}\right]{\theta }^a-\sum _{a=1}^M{\theta }^a=-m,$$

$$G^a:=\sum _{\ell =1}^L{v}^{\ell }{G}^{a,\ell },P_{T_a}:=\sum _{\ell =1}^L{v}^{\ell }{P}_{T_{a,\ell }}$$

Full text

I Introduction

Transfer operator-based methods have attracted lot of attention lately for problems involving dynamical systems analysis and design. In particular, transfer operator-based methods are used for identifying steady state dynamics of the system from the invariant measure of transfer operator, identifying almost invariant sets, and coherent structures [2, 3, 4]. The spectral analysis of transfer operators are also applied for reduced order modeling of dynamical systems with applications to building systems, power grid, and fluid mechanics [5, 6]. Operator-theoretic methods have also been successfully applied to address design problems in control dynamical systems. In particular, transfer operator methods are used for almost everywhere stability verification, control design, nonlinear estimation, and for solving optimal sensor placement problem [7, 1, 8, 9, 10, 11, 12].

In this paper, we continue with the long series of work on the application of transfer operator methods for stability verification and stabilization of nonlinear systems. We develop an analytical and computational framework for the application of transfer operator methods for the stabilization of stochastic nonlinear systems. In [13], we introduced Lyapunov measure for stability verification of stochastic nonlinear systems. We proved that the existence of the Lyapunov measure verifies weaker set-theoretic notion of almost everywhere stochastic stability for discrete-time stochastic systems. Weaker notion of almost everywhere stability was introduced in [14] for continuous time deterministic systems and in [15] for continuous time stochastic systems. In this paper we extend the application of Lyapunov measure for optimal stabilization of stochastic nonlinear systems. Optimal stabilization of stochastic systems is posed as an infinite dimensional linear program. Set-oriented numerical methods are used for the finite dimensional approximation of the transfer operator and the linear program. A key advantage of the proposed transfer operator-based approach for stochastic stability analysis and controller synthesis is that all the stability results along with the computation framework carries over from the deterministic systems [7, 1] to the stochastic systems. The only difference in the stochastic setting is that the transfer Perron-Frobenius operator is defined for the stochastic system.

The results developed in this paper draw parallels from following papers. Lasserre, Hernández-Lerma, and co-workers [16, 17] formulated the control of Markov processes as a solution of the HJB equation. In [18, 19, 20], solutions to stochastic and deterministic optimal control problems are proposed, using a linear programming approach or using a sequence of LMI relaxations. Our paper also draws some connection to research on optimization and stabilization of controlled Markov chains discussed in [21]. Computational techniques based on the viscosity solution of the HJB equation is proposed for the approximation of the value function and optimal controls in [22, Chapter VI].

Our proposed method, in particular the computational approach, draws some similarity with the above discussed references on the approximation of the solution of the HJB equation [23, 24, 22, 25]. Our method, too, relies on discretization of state space to obtain globally optimal stabilizing control. However, our proposed approach differs from the above references in the following two fundamental ways. The first main difference arises due to adoption of non-classical weaker set-theoretic notion of almost everywhere stability for optimal stabilization. The second main difference compared to references [21] and [22] is in the use of the discount factor $\gamma>1$ in the cost function. The discount factor plays an important role in controlling the effect of finite dimensional discretization or the approximation process on the true solution. In particular, by allowing for the discount factor, $\gamma$, to be greater than one, it is possible to ensure that the control obtained using the finite dimensional approximation is truly stabilizing for the nonlinear system [9, 26].

The paper is organized as follows. In section II we present brief overview of results from [13] on Lyapunov measure for stochastic stabilization. In III results on application of Lyapunov measure for optimal stabilization are presented. In section IV, computational framework based on set-oriented numerical methods for finite dimensional approximation of Lyapunov measure and optimal control is presented. Simulation results are presented in section V followed by conclusions in section VI.

II Lyapunov measure for stochastic stability analysis

Consider the discrete-time stochastic system,

[TABLE]

where $x_{n}\in X\subset\mathbb{R}^{d}$ is a compact set. The random vectors, $\xi_{0},\xi_{1},\ldots$, are assumed independent identically distributed (i.i.d) and takes values in $W$ with the following probability distribution,

[TABLE]

and is the same for all $n$ and $v$ is the probability measure. The system mapping $T(x,\xi)$ is assumed continuous in $x$ and for every fixed $x\in X$, it is measurable in $\xi$. The initial condition, $x_{0}$, and the sequence of random vectors, $\xi_{0},\xi_{1},\ldots$, are assumed independent. The basic object of study in our proposed approach to stochastic stability is a linear transfer, the Perron-Frobenius operator, defined as follows:

Definition 1** (Perron-Frobenius (P-F) operator)**

Let ${\cal M}(X)$ be the space of finite measures on $X$. The Perron-Frobenius operator, $\mathbb{P}:{\cal M}(X)\to{\cal M}(X)$, for stochastic dynamical system (1) is given by

[TABLE]

for $\mu\in{\cal M}(X)$, and $A\in{\cal B}(X)$, where ${\cal B}(X)$ is the Borel $\sigma$-algebra on $X$, $T^{-1}_{\xi}(A)=T^{-1}(A,\xi)$ is the inverse image of the set $A$, and $\chi_{A}(x)$ is an indicator function of set $A$.

Assumption 2

We assume $x=0$ is an equilibrium point of system (1), i.e., $T(0,\xi_{n})=0,\;\;\;\forall n,$ for any given sequence of random vectors $\{\xi_{n}\}$.

Assumption 3** (Local Stability)**

We assume the trivial solution, $x=0$, is locally stochastic, asymptotically stable. In particular, we assume there exists a neighborhood ${\cal O}$ of $x=0$, such that for all $x_{0}\in\cal O$,

[TABLE]

and

[TABLE]

where $\xi_{0}^{n}$ notation is used to define the sequence of random variable $\{\xi_{0},\ldots,\xi_{n}\}$.

Assumption 2 is used in the decomposition of the P-F operator in section (II-A) and Assumption 3 is used in the proof of Theorem 7. We will use the notation $U(\epsilon)$ to denote the $\epsilon$ neighborhood of the origin for any positive value of $\epsilon>0$. We have $0\in U(\epsilon)\subset{\cal O}$.

We introduce the following definitions for stability of the stochastic dynamical system (1).

Definition 4** (a.e. stochastic stablity with geometric decay)**

For any given $\epsilon>0$, let $U(\epsilon)$ be the $\epsilon$ neighborhood of the equilibrium point, $x=0$. The equilibrium point, $x=0$, is said to be almost everywhere, almost sure stable with geometric decay with respect to finite measure, $m\in{\cal M}(X)$, if there exists $0<\alpha(\epsilon)<1$, $0<\beta<1$, and $K(\epsilon)<\infty$, such that

[TABLE]

for all sets $B\in{\cal B}(X\setminus U(\epsilon))$, such that $m(B)>0$.

We introduce the following definition of absolutely continuous measures.

Definition 5** (Absolutely continuous measure)**

A measure $\mu$ is absolutely continuous with respect to another measure, $\vartheta$ denoted as $\mu\prec\vartheta$, if $\mu(B)=0$ for all $B\in{\cal B}(X)$ with $\vartheta(B)=0$.

II-A Decomposition of the P-F operator

Let $E=\{0\}$. Hence, $E^{c}=X\setminus E$. We write $T:E\cup E^{c}\times W\to X$. For any set $B\in{\cal B}(E^{c})$, we write

[TABLE]

This is because $T(x,\xi)\in B$ implies $x\notin E$. Since set $E$ is invariant, we define the restriction of the P-F operator on the complement set $E^{c}$. Thus, we can define the restriction of the P-F operator on the measure space ${\cal M}(E^{c})$ as follows:

[TABLE]

for any set $B\in{\cal B}(E^{c})$ and $\mu\in{\cal M}(E^{c})$.

Next, the restriction $T:E\times W\rightarrow E$ can also be used to define a P-F operator denoted by

[TABLE]

where $\mu\in{\cal M}(E)$ and $B\subset{\cal B}(E)$.

The above considerations suggest a representation of the P-F operator, $\mathbb{P}$, in terms of $\mathbb{P}_{0}$ and $\mathbb{P}_{1}$. Indeed, this is the case, if one considers a splitting of the measured space, ${\cal M}(X)={\cal M}_{0}\oplus{\cal M}_{1}$, where ${\cal M}_{0}:={\cal M}(E)$, ${\cal M}_{1}:={\cal M}(E^{c})$, and $\oplus$ stands for the direct sum.

The splitting defined in above equation implies that the P-F operator has a lower-triangular matrix representation given by

[TABLE]

Following definition of Lyapunov measure is introduced for a.e. stability verification of system (1).

Definition 6** (Lyapunov measure)**

A Lyapunov measure, $\bar{\mu}\in{\cal M}(X\setminus U(\epsilon))$, is defined as any positive measure finite outside the $\epsilon$ neighborhood of equilibrium point and satisfies

[TABLE]

for $0<\gamma\leq 1$ and for all sets $B\in{\cal B}(X\setminus U(\epsilon))$.

The following theorem provides the condition for a.e. stochastic stability with geometric decay.

Theorem 7

An attractor set $\cal A$ for the system (1) is a.e. stochastic stable with geometric decay (Definition 4) with respect to finite measure $m$, if and only if for all $\epsilon>0$ there exists a non-negative measure $\bar{\mu}$ which is finite on ${\cal B}(X\setminus U(\epsilon))$ and satisfies

[TABLE]

for all sets $B\subset X\setminus U(\epsilon)$ and for some $\gamma>1$.

Proof:

We omit the proof here due to space constraints. However, the proof follows exactly along the lines of the proof for deterministic systems [1] ∎

III Lyapunov measure for optimal stabilization

We consider the stabilization of stochastic dynamical system

[TABLE]

where $x_{n}\in X\subset\mathbb{R}^{q}$ is the state, $u_{n}\in U\subset\mathbb{R}^{d}$ is the control input, and $\xi_{n}\in W\subset\mathbb{R}^{p}$ is a random variable. The sequence of random variables $\xi_{0},\xi_{1},\ldots$ are assumed to independent identically distributed (i.i.d.) as in (2). For each fixed value of $\xi$ the mapping $T_{\xi}:X\times U\to X$ is assumed to be continuous in $x$ and $u$, and for every fixed values of $x$ and $u$ it is measurable in $\xi$. Both $X$ and $U$ are assumed compact. The objective is to design a deterministic feedback controller, $u_{n}=K(x_{n})$, to optimally stabilize the attractor set $\cal A$.

We define the feedback control mapping $C:X\rightarrow Y:=X\times U$ as $C(x)=(x,K(x))$. We denote by ${\cal B}(Y)$ the Borel-$\sigma$ algebra on $Y$ and ${\cal M}(Y)$ the vector space of real valued measures on ${\cal B}(Y)$. For any $\mu\in{\cal M}(X)$, the control mapping $C$ can be used to define a measure, $\theta\in{\cal M}(Y)$, as follows:

[TABLE]

for all sets $D\in{\cal B}(Y)$ and $B\in{\cal B}(X)$. Since $C$ is an injective function with $\theta$ satisfying (10), it follows from the theorem on disintegration of measure [27] (Theorem 5.8) there exists a unique disintegration $\theta_{x}$ of the measure $\theta$ for $\mu$ almost all $x\in X$, such that $\int_{Y}f(y)d\theta(y)=\int_{X}\int_{C(x)}f(y)d\theta_{x}(y)d\mu(x)$, for any Borel-measurable function $f:Y\to\mathbb{R}$. In particular, for $f(y)=\chi_{D}(y)$, the indicator function for the set $D$, we obtain $\theta(D)=\int_{X}\int_{C(x)}\chi_{D}(y)d\theta_{x}(y)d\mu(x)=[\mathbb{P}_{C}\mu](D).$ Using the definition of the feedback controller mapping $C$, we write

[TABLE]

The system mapping $T:Y\times W\rightarrow X$ can be associated with P-F operators ${\mathbb{P}}_{T}:{\cal M}(Y)\rightarrow{\cal M}(X)$ as

[TABLE]

For the feedback control system $T_{\xi}\circ C:X\times W\rightarrow X$, the P-F operator can be written as a product of ${\mathbb{P}}_{T_{\xi}}$ and ${\mathbb{P}}_{C}$. In particular, we obtain

[TABLE]

[TABLE]

The P-F operators, $\mathbb{P}_{T}$ and $\mathbb{P_{C}}$, are used to define their restriction, $\mathbb{P}_{T}^{1}:{\cal M}({\cal A}^{c}\times U)\to{\cal M}({\cal A}^{c})$, and $\mathbb{P}_{C}^{1}:{\cal M}({\cal A}^{c})\to{\cal M}({\cal A}^{c}\times U)$ to the complement of the attractor set, respectively, similar to Eqs. (5)-(6).

Assumption 8

We assume there exists a feedback controller mapping $C_{0}(x)=(x,K_{0}(x))$, which locally stabilizes the invariant set ${\cal A}$, i.e., there exists a neighborhood $V$ of $\cal A$ such that $T\circ C_{0}(V)\subset V$ and $x_{n}\rightarrow\cal A$ for all $x_{0}\in V$; moreover ${\cal A}\subset U(\epsilon)\subset V$.

Our objective is to construct the optimal stabilizing controller for almost every initial condition starting from $X_{1}$. Let $C_{1}:X_{1}\rightarrow Y$ be the stabilizing control map for $X_{1}(:=X\setminus U(\epsilon))$. The control mapping $C:X\to X\times U$ can be written as follows:

[TABLE]

Furthermore, we assume the feedback control system $T_{\xi}\circ C:X\rightarrow X$ is non-singular with respect to the Lebesgue measure, $m$ for fixed value of $\xi$. We seek to design the controller mapping, $C(x)=(x,K(x))$, such that the attractor set $\cal A$ is a.e. stable with geometric decay rate $\beta<1$, while minimizing the cost function,

[TABLE]

where $x_{0}=x$, the cost function $G:X\times U\times W\to{\mathbb{R}}$ is assumed a continuous non-negative real-valued function for each fixed value of $\xi$ and is assumed to be measurable w.r.t. $\xi$ for fixed values of $x$ and $u$. Furthermore, $G({\cal A},0,\xi)=0$ for all $\xi$, $x_{n+1}=T_{\xi}\circ C(x_{n})$, and $0<\gamma<\frac{1}{\beta}$. The expectation $E_{\xi_{0}^{n}}$ in (12) denotes expectation over the sequence of random variable $\{\xi_{0},\xi_{1},\ldots,\xi_{n}\}$. Note, that in the cost function (12), $\gamma$ is allowed to be greater than one and this is one of the main departures from the conventional optimal control problem, where $\gamma\leq 1$. Under the assumption that the controller mapping $C$ renders the attractor set a.e. stable with a geometric decay rate, $\beta<\frac{1}{\gamma}$, the cost function (12) is finite.

Remark 9

We will use the notion of the scalar product between continuous function $h\in{\cal C}^{0}(X)$ and measure $\mu\in{\cal M}(X)$ as $\left<h,d\mu\right>_{X}:=\int_{X}h(x)d\mu(x)$ [28].

Let the controller mapping, $C(x)=(x,K(x))$, be such that the attractor set $\cal A$ for the feedback control system $T_{\xi}\circ C:X\to X$ is a.e. stable with geometric decay rate $\beta<1$. Then, the cost function (12) is well defined for $\gamma<\frac{1}{\beta}$ and, furthermore, the cost of stabilization of the attractor set $\cal A$ with respect to Lebesgue almost every initial condition starting from the set $B\in{\cal B}(X_{1})$ can be expressed as follows:

[TABLE]

where, ${\cal Z}={\cal A}^{c}\times U\times W$, $x_{0}=x$ and $\bar{\mu}_{B}$ is the solution of the following control Lyapunov measure equation,

[TABLE]

for all $D\in{\cal B}(X_{1})$ and where $m_{B}(\cdot):=m(B\cap\cdot)$ is a finite measure supported on the set $B\in{\cal B}(X_{1})$.

The minimum cost of stabilization is defined as the minimum over all a.e. stabilizing controller mappings, $C$, with a geometric decay as follows:

[TABLE]

Next, we write the infinite dimensional linear program for the optimal stabilization of the attractor set $\cal A$. Towards this goal, we first define the projection map, $P_{1}:{\cal A}^{c}\times U\rightarrow{\cal A}^{c}$ as: $P_{1}(x,u)=x,$ and denote the P-F operator corresponding to $P_{1}$ as $\mathbb{P}_{P_{1}}:{\cal M}({\cal A}^{c}\times U)\rightarrow{\cal M}({\cal A}^{c})$, which can be written as $[{\mathbb{P}}^{1}_{P_{1}}\theta](D)=\int_{{\cal A}^{c}\times U}\chi_{D}(P_{1}(y))d\theta(y)=\int_{D\times U}d\theta(y)=\mu(D)$. Using this definition of projection mapping, $P_{1}$, and the corresponding P-F operator, we can write the linear program for the optimal stabilization of set $B$ with unknown variable $\theta$ as follows:

[TABLE]

for $D\in{\cal B}(X_{1})$.

Remark 10

Observe the geometric decay parameter satisfies $\gamma>1$. This is in contrast to most optimization problems studied in the context of Markov-controlled processes, such as in Lasserre and Hernández-Lerma [16]. Average cost and discounted cost optimality problems are considered in [16, 22]. The additional flexibility provided by $\gamma>1$ guarantees the controller obtained from the finite dimensional approximation of the infinite dimensional program (16) also stabilizes the attractor set for control dynamical system.

IV Computational approach

We discretize the state-space and control space for the purposes of computations as described below. Borrowing the notation from [9], let ${\cal X}_{N}:=\{D_{1},...,D_{i},...,D_{N}\}$ denote a finite partition of the state-space $X\subset\mathbb{R}^{q}$. The measure space associated with ${\cal X}_{N}$ is $\mathbb{R}^{N}$. We assume without loss of generality that the attractor set, ${\cal A}$, is contained in $D_{N}$, that is, ${\cal A}\subseteq D_{N}$. The control space, $U$, is quantized and the control input is assumed to take only finitely many control values from the quantized set, ${\cal U}_{M}=\{u^{1},\ldots,u^{a},\ldots,u^{M}\}$, where $u^{a}\in\mathbb{R}^{d}$. The partition, ${\cal U}_{M}$, is identified with the vector space, $\mathbb{R}^{d\times M}$. Similarly, the space of uncertainty, $W$, and the probability measure $v$ is quantized and are assumed to take only finitely many values ${\cal W}=\{\xi^{1},\ldots,\xi^{\ell},\ldots,\xi^{L}\}$, and $\vartheta=\{v^{1},\ldots,v^{\ell},\ldots,v^{L}\}$ where $\xi^{\ell}\in\mathbb{R}^{p}$ and $0\leq v^{\ell}\leq 1$ for all $\ell$ and $\sum_{\ell=1}^{L}v^{\ell}=1$. The discrete probability measure on the finite dimensional uncertainty space is assigned as follows:

[TABLE]

The space of uncertainty is identified with finite dimensional space $\mathbb{R}^{p\times L}$. The system map that results from choosing the controls $u=u^{a}$ and uncertainty value $\xi=\xi^{\ell}$ is denoted by $T_{u^{a},\xi^{\ell}}$ and the corresponding P-F operator is denoted as $P_{T_{u^{a},\xi^{\ell}}}\in\mathbb{R}^{N\times N}$. Note that for system mapping $T_{u^{a},\xi^{\ell}}$, the control on all sets of the partition is $u(D_{i})=u^{a}$, for all $D_{i}\in{\cal X}_{N}$. For brevity of notation, we will denote the P-F matrix $P_{T_{u^{a},\xi^{\ell}}}$ by $P_{T_{a,\ell}}$ and its entries are calculated as

[TABLE]

where $m$ is the Lebesgue measure and $(P_{T_{a,\ell}})_{(ij)}$ denotes the $(i,j)$-th entry of the matrix. Since $T_{u^{a},\xi^{\ell}}:X\rightarrow X$, we have $P_{T_{a,\ell}}$ is a Markov matrix. Additionally, $P^{1}_{T_{a,\ell}}:\mathbb{R}^{N-1}\rightarrow\mathbb{R}^{N-1}$ will denote the finite dimensional counterpart of the P-F operator restricted to ${\cal X}_{N}\setminus D_{N}$, the complement of the attractor set. It is easily seen that $P^{1}_{T_{a,\ell}}$ consists of the first $(N-1)$ rows and columns of $P_{T_{a,\ell}}$.

With the above quantization of the control space and partition of the state space, the determination of the control $u(x)\in U$ (or equivalently $K(x)$) for all $x\in{\cal A}^{c}$ has now been cast as a problem of choosing $u_{N}(D_{i})\in{\cal U}_{M}$ for all sets $D_{i}\subset{\cal X}_{N}$. The finite dimensional approximation of the optimal stabilization problem (16) is equivalent to solving the following finite-dimensional LP:

[TABLE]

where we have used the notation $(\cdot)^{{}^{\prime}}$ for the transpose operation, $m\in\mathbb{R}^{N-1}$ and $(m)_{(j)}>0$ denote the support of Lebesgue measure, $m$, on the set $D_{j}$, $G^{a,\ell}\in\mathbb{R}^{N-1}$ is the cost defined on ${\cal X}_{N}\setminus D_{N}$ with $(G^{a,\ell})_{(j)}$ the cost associated with using control action $u^{a}$ on set $D_{j}$ with uncertainty value $\xi=\xi^{\ell}$; $\theta^{a}\in\mathbb{R}^{N-1}$ are, respectively, the discrete counter-parts of infinite-dimensional measure quantities in (16). We define following quantities

[TABLE]

to rewrite finite-dimensional LP (17) as follows:

[TABLE]

In the LP (18), we have not enforced the constraint,

[TABLE]

for each $j=1,...,(N-1)$. The above constraint ensures the control on each set in unique. We prove in the following the uniqueness can be ensured without enforcing the constraint, provided the LP (18) has a solution. To this end, we introduce the dual LP associated with the LP in (18). The dual to the LP in (18) is,

[TABLE]

In the above LP (20), $V$ is the dual variable to the equality constraints in (18).

IV-A Existence of solutions to the finite LP

We make the following assumption throughout this section.

Assumption 11

There exists $\theta^{a}\in\mathbb{R}^{N-1}\;\forall\;a=1,\ldots,M$, such that the LP in (18) is feasible for some $\gamma>1$.

Lemma 12

Consider a partition ${\cal X}_{N}=\{D_{1},\ldots,D_{N}\}$ of the state-space $X$ with attractor set ${\cal A}\subseteq D_{N}$ and a quantization ${\cal U}_{M}=\{u^{1},\ldots,u^{M}\}$ of the control space $U$. Suppose Assumption 11 holds for some $\gamma>1$ and for $m,G>0$. Then, there exists an optimal solution, $\theta$, to the LP (18) and an optimal solution, $V$, to the dual LP (20) with equal objective values, ($\sum\limits_{a=1}^{M}(G^{a})^{{}^{\prime}}\theta^{a}=m^{\prime}V$) and $\theta,V$ bounded.

Proof:

Refer to proof of Lemma 12 in [1]. ∎

The next result shows the LP (18) always admits an optimal solution satisfying (19).

Lemma 13

Given a partition ${\cal X}_{N}=\{D_{1},\ldots,D_{N}\}$ of the state-space, $X$, with attractor set, ${\cal A}\subseteq D_{N}$, and a quantization, ${\cal U}_{M}=\{u^{1},\ldots,u^{M}\}$, of the control space, $U$. Suppose Assumption 11 holds for some $\gamma>1$ and for $m,G>0$. Then, there exists a solution $\theta\in\mathbb{R}^{N-1}$ solving (18) and $V\in\mathbb{R}^{N-1}$ solving (20) for any $\gamma\in[1,\overline{\gamma}_{N})$. Further, the following hold at the solution: 1) For each $j=1,...,(N-1)$, there exists at least one $a_{j}\in 1,...,M$, such that $(V)_{(j)}=\gamma(P^{1}_{T_{a_{j}}}V)_{(j)}+(G^{a_{j}})_{(j)}$ and $(\theta^{a_{j}})_{(j)}>0$. 2) There exists a $\tilde{\theta}$ that solves (18), such that for each $j=1,...,(N-1)$, there is exactly one $a_{j}\in 1,...,M$, such that $(\tilde{\theta}^{a_{j}})_{(j)}>0$ and $(\tilde{\theta}^{a^{\prime}})_{(j)}=0$ for $a^{\prime}\neq a_{j}$.

Proof:

Refer to proof of Lemma 14 in [1]. ∎

The following theorem states the main result.

Theorem 14

Consider a partition ${\cal X}_{N}=\{D_{1},\ldots,D_{N}\}$ of the state-space, $X$, with attractor set, ${\cal A}\subseteq D_{N}$, and a quantization, ${\cal U}_{M}=\{u^{1},\ldots,u^{M}\}$, of the control space, $U$. Suppose Assumption 11 holds for some $\gamma>1$ and for $m,G>0$. Then, the following statements hold: 1) there exists a bounded $\theta$, a solution to (18) and a bounded $V$, a solution to (20); 2) the optimal control for each set, $j=1,...,(N-1)$, is given by $u(D_{j})=u^{a(j)},\;\;\mbox{where }a(j):=\min\{a|(\theta^{a})_{(j)}>0\}$; 3) $\mu$ satisfying $\gamma(P^{1}_{T_{u}})^{{}^{\prime}}\mu-\mu=-m\;\;,\mbox{ where }(P^{1}_{T_{u}})_{(ji)}=(P^{1}_{T_{a(j)}})_{(ji)}$ is the Lyapunov measure for the controlled system.

Proof:

Refer to proof of Lemma 15 in [1]. ∎

V Example: Inverted Pendulum on cart

[TABLE]

where $g=9.8,l=0.5,m=2,M=8,\zeta=0,m_{r}=\cfrac{m}{(m+M)},a=\cfrac{g}{l},b=\cfrac{m_{r}}{ml}$. The cost function is assumed to be $G(x,u)=x^{2}+\dot{x}^{2}+u^{2}$. For uncontrolled system, $u=0$,

there are two equilibrium points, one equilibrium point at $(\pi,0)$ is stable in Lyapunov sense with eigenvalues of linearization on the $j\omega$ axis, the second equilibrium point at the origin is a saddle and unstable. In Fig. 1 we show the phase portrait for the uncontrolled system. The objective is to optimally stabilize the saddle equilibrium point at the origin. For the purpose of discretization we use $\delta t=0.1$ as time discretization for the simulations. The state space $X$ is chosen to be limited in $[-\pi,\pi]\times[-10,10]$ and is partitioned into $70\times 70=4900$ boxes. For constructing the P-F matrix $10$ initial conditions are located in each box. The control set is discretized as follows ${\cal U}=\{-80,-70\ldots,-10,0,10,\ldots,70,80\}$.

• •

Case 1: The damping parameter $\zeta$ is assumed to be random and uniformly distributed with mean zero and uniformly supported on the interval $[-\sigma,\sigma]$. Similarly, the range of random parameter $[-\sigma,\sigma]$ is divided into 10 uniformly spaced discrete values for random parameter $\xi$.

• •

Case 2: The parameter $b$ multiplying the control input is assumed to be Bernoulli random variable with statistics $Prob(b=1)=p$ and $Prob(b=0)=1-p$ for every time.

In Fig. 2, we show the plot for the Lyapunov measure, optimal control for $\sigma=0.1$. In Fig. 3, we show the plot for the optimal cost and the percentage of initial condition that are attracted to the origin. It is interesting to notice that the optimal cost along the stable manifold of the uncontrolled system is small whereas along the unstable manifold is large. This is because of the fact the optimal control is design to exploit the natural dynamics of the system since there is non-zero cost on the control efforts. The simulation result in Fig.3b are obtained by performing time domain simulation with eight initial conditions in each box iterated over $100$ time step. We notice that close to $100$ percentage of initial conditions are attracted to the origin. In Fig. 4 we show the comparison of sample trajectory for the open loop and closed loop system. The sample trajectory shows that feedback controller is able to stabilize the origin.

Case 2: For case 2, we consider Bernoulli uncertainty for the input channel. In this case the random parameter $b$ can take only two values at every time instant i.e., $b=0$ or $b=1$. For this case we only show the plots for the percentage of initial conditions that can be optimally stabilized to the origin for two different values of erasure probabilities $1-p=0.15$ and $1-p=0.5$. Again the simulation results for this case are obtained by performing time domain simulation with eight initial conditions in each box iterated over $100$ time step. From Fig. 5, we notice that with erasure probability of $1-p=0.15$ more than $97\%$ of initial conditions are attracted to the origin. While for $1-p=0.5$ only $66\%$ of points are attracted to the origin thereby indicating that the origin is not stabilized with erasure probability of $0.5$.

VI Conclusions

Transfer Perron-Frobenius operator-based framework is introduced for optimal stabilization of stochastic nonlinear systems. Weaker set-theoretic notion of almost everywhere stability is used for the design of optimal stabilizing feedback controller. The optimal stabilization problem is formulated as an infinite dimensional linear program. The finite dimensional approximation of the linear program and the associated optimal feedback controller is obtained using set-oriented numerics.