# A Stochastic Optimal Regulator for a Class of Nonlinear Systems

Authors: Gabriella Mavelli; Giovanni Palombo; Pasquale Palumbo
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

This work investigates an optimal control problem for a class of stochastic differential bilinear systems, affected by a persistent disturbance provided by a nonlinear stochastic exogenous system (nonlinear drift and multiplicative state noise). The optimal control problem aims at minimizing the average value of a standard quadratic-cost functional on a finite horizon. It has been supposed that neither the state of the system nor the state of the exosystem is directly measurable (incomplete information case). The approach is based on the Carleman embedding, which allows to approximate the nonlinear stochastic exosystem in the form of a bilinear system (linear drift and multiplicative noise) with respect to an extended state that includes the state Kronecker powers up to a chosen degree. This way the stochastic optimal control problem may be restated in a bilinear setting and the optimal solution is provided among all the affine transformations of the measurements. The present work is a nontrivial extension of previous work of the authors, where the Carleman approach was exploited in a framework where only additive noises had been conceived for the state and for the exosystem. Numerical simulations support theoretical results by showing the improvements in the regulator performances by increasing the order of the approximation.


---

## Body

## 1. Introduction

Consider the following stochastic differential system described by means of the Itô formalism:(1)dxt=Axtdt+Butdt+Lztdt+∑i=1bMi+NixtdWixt,x0=x0,dyt=Cxtdt+GdWyt,y0=0a.s.,where xt∈Rn is the state of the system, ut∈Rp is the control input, and yt∈Rq is the measured output. Wxt=W1x…WbxT∈Rb and Wyt∈Rd are independent standard Wiener processes with respect to a family of increasing σ-algebras ℱt,t≥0, referred to a probability space Ω,ℱ,P. The standard assumption rank GGT=q is made. The initial state x0 is an ℱ0-measurable random vector, independent of Wxt and Wyt. zt∈Rm is a persistent disturbance generated by the following nonlinear stochastic exogenous differential system (the exosystem):(2)dzt=ϕztdt+∑j=1hFj+DjztdWjzt,z0=z0,where ϕ:Rm⟼Rm is a smooth nonlinear map. Wjzt is the j‐th element of Wzt∈Rh which represents a standard Wiener process with respect to ℱt,t≥0, independent of the state and output noise processes Wxt and Wyt.Linear and nonlinear exosystems have been widely exploited to model uncertainties as well as sustained perturbations, especially within applicative engineering frameworks such as missile systems, robotics, and wind turbines (one can refer to [1] and references therein, where an observer-based approach is exploited to estimate the unknown exosystem). Explicit noises in the exosystem dynamics could further enhance the nonlinearities of the system, and they have been considered in the recent literature for an exosystem with linear drift [1].The optimal control problem here investigated refers to the following standard quadratic-cost index to be minimized on a finite horizon0,tf:(3)Ju=12ExTtfSxtf+∫0tfxTtQxt+uTtRutdt,where S, Q are symmetric positive-semidefinite matrices and R a symmetric positive-definite matrix. E⋅ stands for the expectation value operator.The problem under investigation is clearly framed in the context of stochastic optimal control problems for nonlinear systems. Even though neglecting the stochastic disturbances, the nonlinear fashion of the optimal control problem does not ensure analytical solutions from the application of the maximum principle, because it requires the solution of the Two-Point Boundary Value problem (see, e.g., [2–4] and references therein). Dealing with stochastic systems, such investigation has been usually carried out in the framework of a complete knowledge of the state of the system (i.e., the state is directly available with no need for outputs providing possibly noisy/incomplete measurements of the state). Usual difficulties involve the solution of the Hamilton–Jacobi–Bellman (HJB) equations associated to the optimal control problem: in [5], the stochastic HJB equation is iteratively solved with successive approximations; in [6], the infinite-time HJB equation is reformulated as an eigenvalue problem; in [7], a transformation approach is proposed for solving the HJB equation arising in quadratic-cost control for nonlinear deterministic and stochastic systems. Finally, in a pair of recent papers, a solution to the nonlinear HJB equation is provided, by expressing it in the form of decoupled Forward and Backward Stochastic Differential Equations (FBSDEs), for an L2- and an L1-type optimal control setting (see [8, 9], respectively). As stated above, the solutions proposed in these references rely on a complete knowledge of the state of the system; thus, they do not require any nonlinear state-estimation algorithm to infer information from noisy measurements. To the best of our knowledge, the only reference that deals with stochastic optimal control problems in a nonlinear framework with incomplete knowledge of the state is [10], though nonlinearities are restricted only to the diffusion term where the noise affects the state dynamics in a multiplicative fashion; the state drift and the output equation providing noisy measurements are both linear.To cope with the incomplete information case, a state-estimation algorithm is required. The optimal state estimate among all the Borel transformations of the measurements, in this case, requires the knowledge of the whole conditional probability density provided by the solution of the Kolmogorov equation, a nontrivial infinite-dimensional problem. Several methods can be found in the literature in order to achieve it, dealing with techniques inherently based on the searching for PDE numerical solutions (see, e.g., the recent approaches on finite elements methods [11, 12]) or with Monte Carlo approaches such as, among the others, particle methods [13] or multilevel Monte Carlo methods [14]. All these approaches share a nontrivial computational cost.A different philosophy consists in introducing an approximation of the original setting according to which the optimization problem is restated in a form for which there exist available solutions in the literature. In this case, a tradeoff should be searched between the simplifications provided by the approximation and its displacement from the real case. For instance, the Extended Kalman Filter relies on the linearization of a stochastic nonlinear system and is among the most widely used algorithms for real-time state estimate because of its simplicity [15]; nonetheless, there are many applications where linearization is a very coarse approximation of reality and filters simply do not work.In the spirit of the aforementioned philosophy, in this note, we apply the Carleman approximation to the nonlinear exosystem. The Carleman approach consists in the embedding of the original nonlinear differential stochastic system onto an infinite-dimensional system whose state accounts for the Kronecker powers of any order of the original state. With respect to such a state, the dynamics can be written in a bilinear fashion (linear drift and multiplicative noise), and theν-degree Carleman approximation is achieved by truncating the higher-than-ν Kronecker powers. The idea is further supported by the new results on polynomial filtering, which take advantage of the polynomial structure of the problem to achieve more accurate estimations, as in [16, 17]. Bilinear systems have gained an increasing interest since early seventies, when they have been started to be investigated as an appealing class of nearly linear systems, [18]; according to more recent literature, they can be found in different fields of engineering and mathematical sciences, including economics, electronic circuits, and theoretical biology (see [19] and references therein). Moreover, there can be found suboptimal state-estimation algorithms, suitably designed for stochastic bilinear systems [20, 21]. Within this framework, the Carleman embedding technique has been successfully applied in the recent years both to a discrete- and continuous-time framework to solve filtering problems by first reformulating them in a bilinear fashion and, then, by applying known suboptimal algorithms (see, e.g. [22–25]).In this note, once the Carleman bilinear approximation of the exosystem is coupled to the state equations, the optimal solution of the reformulated problem is still not implementable because of the nonavailability of a finite-dimensional algorithm for the optimal control. Therefore, we propose theoptimal linear regulator, by extending the results of [10], consisting of the optimal solution among all of the Rp-valued square-integrable affine transformations of the observations. A similar approach can also be found in [26] in a discrete-time framework, and in [27] in a continuous-time framework, where only additive noises had been conceived for the state and for the exosystem. The present note is actually a nontrivial extension of [27], since multiplicative noises are now considered both in the system dynamics and in the exosystem one.


## 2. Carleman Approximation of the Stochastic Exosystem

Consider the Taylor polynomial expansion around a given pointz˜∈Rm for the exosystem, supposed to exist according to standard analyticity hypotheses. By properly exploiting the Kronecker formalism, defining the displacement φt=zt−z˜, it follows:(4)ϕzt=∑i=0∞Aiz˜φti,(5)Aiz˜=∇zi⊗ϕzti!zt=z˜.The square brackets denote the Kronecker power, and the differential operator∇zi⊗ applied to a generic function ψ=ψz:Rm⟼Rm is defined as follows:(6)∇z0⊗ψ=ψ,∇zi+1⊗ψ=∇z⊗∇zi⊗ψ,i≥0,with ∇z=∂/∂z1…∂/∂zm and ∇z⊗ψ the Jacobian of the vector function ψ (see the Appendix for a quick survey on the Kronecker algebra). Thus, by taking into account (2) and (4), we have that(7)dφt=dzt=∑i=0∞Aiz˜φtidt+∑j=1hFj+Djz˜+DjφtdWjzt=∑i=0∞Aiz˜φtidt+∑j=1hF¯jz˜+DjφtdWjzt,where F¯jz˜=Fj+Djz˜. We will drop hereafter the explicit dependence of z˜ in F¯jz˜ to shorten notation. The differential of the Kronecker powers of the displacement in (7) is then required to be computed in order to build up the Carleman embedding. By standard Itô calculus [28], it follows, for any k>1,(8)dφtk=∇φ⊗φtk∑i=0∞Aiz˜φtidt+12∇φ2⊗φtk∑j=1hF¯jz˜+Djφt2dt+∇φ⊗φtk∑j=1hF¯jz˜+DjφtdWjzt.By exploiting Lemma 1 in the Appendix, according to the definition of matricesUmk and Omk in (A.4),(9)dφtk=∑i=k−1∞A¯k,iφtidt+12∑j=1hOmkF¯j2⊗Imk−2φtk−2+C¯j⊗Imk−2φtk−1+Dj2⊗Imk−2φtkdt+∑j=1hUmkF¯j⊗Imk−2φtk−1+Dj⊗Imk−2φtkdWjzt,with Im denoting the identity matrix in Rm×m and(10)A¯k,i=UmkAi−k+1z˜⊗Imk−1,C¯j=F¯j⊗Dj+Dj⊗F¯j.By means of the definitions(11)ℱj,k=UmkF¯j⊗Imk−1,ℬj,k=UmkDj⊗Imk−1,F˜h,k=12∑j=1hOmkF¯j2⊗Imk−2,C˜h,k=12∑j=1hOmkC¯j⊗Imk−2,B˜h,k=12∑j=1hOmkDj2⊗Imk−2,equation (8) can be written in more compact form:(12)dφtk=∑i=k−1∞A¯k,iz˜φti+F˜h,kφtk−2+C˜h,kφtk−1+B˜h,kφtkdt+∑j=1hℱj,kφtk−1+ℬj,kφtkdWjzt.Collecting the Kronecker powersφtk, k=1,2,… in a unique vector Φt, one obtains the following infinite-dimensional bilinear differential equation:(13)dΦt=Λz˜Φt+Θz˜dt+∑j=1hΓjΦt+Γj,0dWjzt,with Φt=φtTφt2T…T and(14)Λz˜=A¯1,1z˜A¯1,2z˜A¯1,3z˜A¯1,4z˜⋯ϒ2z˜Ω2z˜A¯2,3z˜A¯2,4z˜⋯F˜h,3ϒ3z˜Ω3z˜A¯3,4z˜⋯0m4×mF˜h,4ϒ4z˜Ω4z˜⋯⋮⋮⋮⋮⋮⋱,with(15)ϒ2z˜=A¯2,1z˜+C˜h,2,Ω2z˜=A¯2,2z˜+B˜h,2,ϒ3z˜=A¯3,2z˜+C˜h,3,Ω3z˜=A¯3,3z˜+B˜h,3,ϒ4z˜=A¯4,3z˜+C˜h,4,Ω4z˜=A¯4,4z˜+B˜h,4,(16)Θz˜=A¯1,0Tz˜F˜h,2T01×m301×m4…T,(17)Γj=ℬj,10m×m20m×m30m×m4⋯ℱj,2ℬj,20m2×m30m2×m4⋯0m3×mℱj,3ℬj,30m3×m4⋯0m4×m0m4×m2ℱj,4ℬj,4⋯⋮⋮⋮⋮⋱,(18)Γj,0=ℱj,1T01×m201×m301×m4…T.Theν-degree Carleman approximation consists in collecting the first ν components of vector Φt in the finite-dimensional vector:(19)Φνt=φtTφt2T…φtνTT,and, then, describing its dynamics according to the finite-dimensional version of (13):(20)dΦνt=Λ˜z˜Φνt+Θ˜z˜dt+∑j=1hΓ˜jΦνt+Γ˜j,0dWjzt,in which Λ˜z˜, Θ˜z˜, Γ˜j, Γ˜j,0 are the finite-dimensional matrices achieved by accounting for the first ν blocks of the Carleman embedding matrices, i.e., accounting for the first mν=m+m2+⋯+mν rows and columns of (14)–(18). Then, we need to substitute zt in (1) by means of its Carleman approximation provided by Φνt. By doing this, the state dynamics no more refers to the original xt; therefore, it will be replaced by xνt. For the same reason, we replace yt with yνt.Remark 1.
It is worth noting that the Itô correction term in (8) introduces nontrivial blocks in Λ, the dynamic matrix of the Carleman embedding, defined in (14)-(15), see also (11). These blocks could play an active role in determining the stability properties of the Carleman approximation. Clearly, such investigation gains much more importance for control problems that involve the asymptotic behavior of the system, such as in optimal control with infinite-horizon cost functionals.
Therefore, the state dynamics (1), endowed with the output measurements, are now replaced by the following equations:(21)dxνt=Axνtdt+Butdt+L˜Φνt+z˜νdt+∑i=1bMi+NixνtdWixt,dyνt=Cxνtdt+GdWyt,where(22)L˜=L0n×m2+⋯+mν,z˜ν=z˜T01×m2+⋯+mνT.
Finally, with the aim of writing equation (22) in a more compact form, we define the extended state vector(23)χνt=xνtΦνt,according to which the finite-dimensional bilinear extended system is achieved:(24)dχνt=Aνz˜χνt+ℬνut+Dνz˜dt+∑j=1b+hℋν,jχνt+ℱν,jdWt,dyνt=Cνχνt+GdWyt,(25)Aνz˜=AL˜0mν×nΛ˜z˜,ℬν=B0mν×p,(26)ℋν,j=Nj0n×mν0mν×n0mν×mν,j=1,2,…,b,0n×n0n×mν0mν×nΓ˜j−b,j=b+1,b+2,…,b+h,(27)ℱν,j=Mj0mν×1,j=1,2,…,b,0n×1Γ˜j−b,0,j=b+1,b+2,…,b+h,(28)Dνz˜=L˜z˜νΘ˜z˜,Cν=C0q×mν,Wt=WtxWtz.
Considering (24), the cost functional (3) becomes(29)Jνu=12EχνTtfSνχνtf+∫0tfχνTtQνχνt+uTtRutdt,with(30)Sν=S0n×mν0mν×n0mν×mν,Qν=Q0n×mν0mν×n0mν×mν.


## 3. Optimal Linear Regulator

By means of the Carleman approximation scheme, the original nonlinear optimal control problem (1)–(3) is now restated in the problem of minimizing (29), subject to the bilinear system (24). As already stated in the Introduction, the optimal solution is still not affordable according to a finite-dimensional state-estimate algorithm; therefore, we look for suboptimal solutions. To this end, we synthesize the solution ut providing the minimum of index (29) among all the Rp-valued square-integrable affine transformations of the random variables yντ,0≤τ≤t≤tf. Such a problem has been properly formalized (quadratic functional cost and bilinear differential system) in [10], where a solution is given. It is worth noticing that the solution provided in [10] is not straightforwardly applicable here because of a constant deterministic drift and of an additive noise in the state equation of the extended system. The extension of [10] to such a case has been presented in [27]: indeed, although the original nonlinear frameworks addressed in this note and in [27] are different, the Carleman approximation provides the same mathematical structure for both the embedded system. This fact highlights the advantages of the Carleman approximation, that allows to restate quite different filtering/control problems into the unifying bilinear formulation. Besides, in the spirit of [10], the results proposed can be interpreted as a Separation Principle in a suboptimal sense, since the optimal linear filter is designed independently of the optimal regulator that benefits of the state estimate according to the incomplete information case. In summary, results in [27] are quite straightforwardly applicable and resumed in the following theorem that somehow resembles the Separation Principle.Theorem 1 [27].
Suppose a solution exists for the following backward generalized Riccati equations:(31)V˙t=−AνTz˜Vt−VtAνz˜−Qν−∑k=1b+hℋνkTVtℋνk+VtℬνR−1ℬνTVt,Vtf=Sν,g˙t=−2VtDνz˜−AνTz˜gt+VtℬνR−1ℬνTgt−2∑k=1b+hℋνkTVtℱνk,gtf=0,with Vt=VTt≥0. Then, the solution to the optimal control problem of minimizing the cost criterion (29), under the differential constraints (24), with ut∈ℒtpyν is given by(32)uot=Lotχν^t+αt,with(33)Lot=−R−1ℬνTVt,αt=−12R−1ℬνTgt,where χν^t is the optimal (in the sense of the minimum error variance) estimate of χνt among all the Rn+mν-valued square-integrable affine transformations of yντ,0≤τ≤t≤tf, which is the projection of χν onto ℒtn+mνyν:(34)χν^t=Πχνtℒtn+mνyν,(formally the projection onto ℒtn+mνyν is a random variable such that the difference χνt−Πχνtℒtn+mνyν is orthogonal to ℒtn+mνyν, i.e., is uncorrelated with all random variables in ℒtn+mνyν).Remark 2.
It is worth noting that sums in equation (31) include more terms than the corresponding ones exploited in [27]. This is because, different from [27], here, the multiplicative noise in the state dynamics makes nontrivial matrices ℋν,j for j=1,…,b, see (26).
Concerning the optimal linear filter providingχν^t in (32), the following theorem provides the equations. Its proof is a straightforward consequence of Theorem 4.2 in [10].Theorem 2.
Consider the stochastic system (24) with ut=uot=Lotχν^t+αt, Lot, αt as in (33), and χν^t as in (34). Then, χν^t satisfies the equation:(35)dχν^t=Aνz˜+ℬνLotχν^t+ℬναt+Dνz˜dt+Ktdyνt−Cνχν^tdt,with χν^0=Eχν0 and Kt=PtCνTGGT−1, where Pt is the error covariance matrix evolving according to the equation:(36)P˙t=Aνz˜Pt+PtAνTz˜+∑j=1b+hℋν,jΨνtℋν,jT+∑j=1b+hℋν,jμνt+ℱν,jℋν,jμνt+ℱν,jT−PtCνTGGT−1CνPt,with P0=Covχν0 and(37)μνt=Eχνt,Ψνt=Covχνt,obeying the following equations:(38)μ˙νt=Aνz˜+ℬνLotμνt+ℬναt+Dνz˜,Ψ˙νt=Aνz˜+ℬνLotΨνt+ΨνtAνz˜+ℬνLotT+∑j=1b+hℋν,jΨνtℋν,jT+∑j=1b+hℋν,jμνt+ℱν,jℋν,jμνt+ℱν,jT.Note that according to the optimal initialization of the filtering algorithm associated to the proposed control law, the second-order moments of the initial extended stateχν0 have to be finite and available, which means finite and available moments up to order 2 for x0 and up to order 2ν for z0.Remark 3.
The filter proposed in Theorem 2 provides the optimal linear estimate ofχν as a function of the observations yν. However, the available measurements are given by the output y (instead of yν); therefore, the differential dyν in (33) should be replaced by dy.


## 4. Numerical Simulations

Numerical simulations refer to a second-order system for (1) with scalar input and scalar output:(39)A=−4.5−20.6−1.5,B=12,L=12−13,M1=0.122,N1=0.1011−2,C=10,M2=0.11−1,N2=0.1−0.510.22,G=2.According to (2), two kinds of exosystems are considered:(40)ϕ1z=−z1+z1z2z1−2z2−z1z2,ϕ2z=z2+0.1z1z2−z1.Disregarding the noise,ϕ1 provides a unique asymptotically stable equilibrium in the origin. The linear approximation around the origin exhibits the same qualitative behavior of the nonlinear exosystem, since it is asymptotically stable. As a matter of fact, the first-order Carleman approximation is expected to work well, with higher-order Carleman approximations playing a marginal role. On the other hand, with regard to ϕ2, we have a unique unstable equilibrium point in the origin, and the qualitative behavior exhibits a limit cycle in absence of noise. By applying the Carleman linearization around the origin, we find that the approximate linear exosystem is marginally stable. As a matter of fact, for both the original nonlinear system and for its linear approximation, disregarding the noises, we have sustained oscillations converging to different limit cycles (the one related to the linear approximation that strongly depends of the initial conditions). This fact may result in an unsatisfying linear approximation, and the addition of noises may further worsen the performances.The following matrices are involved in the exosystem (2):(41)F1=0.11−2,D1=0.1−120−2,F2=0.1−11,D2=0.1−0.5024.Initial conditions arex0=12T, z0=0.40T for all the simulations, whilst initial conditions for the filtering algorithm are x^0=34T, z^0=3.43T. The cost functional weight matrices are S=2⋅I2×2, Q=5⋅I2×2, and R=2. The Carleman approximation of the exosystem is achieved around the origin.100 numerical simulations have been produced for both the types of exosystems, by using the Euler–Maruyama method [29], with step integration δ=10−3 in the time interval 0,5. As expected, with regard to ϕ1, results show increasing improvements of the performance index Jν up to degree ν=4, although according to little improvements with respect to the first-order case (J1=5.28±0.02, J2=4.96±0.02, J3=4.91±0.02, J4=4.90±0.02). Results are shown in Figure 1. On the other hand, with regards to ϕ2, results show that the second-order Carleman approximation is enough to provide a significant improvement in the cost functional (J1=202.1±3.3, J2=173.1±3.6) with a reduction of the average value of the cost functional of about 15%. Higher-order approximations exhibit a behavior analogous to the second-order case (the cost functional seems to approach a plateau). Results are shown in Figure 2.Figure 1
Values ofJν for ν=1,2,3,4 for the system with ϕ1z.Figure 2
Values ofJν for ν=1,2,3,4 for the system with ϕ2z.


## 5. Conclusions

A stochastic optimal control problem has been investigated for bilinear stochastic differential systems, driven by a persistent perturbation provided by an exogenous stochastic nonlinear system with multiplicative noise. This work represents a nontrivial extension of a previous work where additive noise (instead of multiplicative noise as in the present case) was considered affecting both state and exosystem. The approach followed relies on the Carleman embedding, successfully exploited in the stochastic framework in the last decade both for filtering and control purposes.



---

*Source: 9763193-2019-10-10.xml*