Second variation of trajectory of a controlled dynamical system

An integral representation for the second variation of trajectory of a dynamical system under control is obtained. This representation contains some tensor of the third rank introduced here. A differential equation for this tensor is presented. A second order method for solution of the optimal control problem based on the second variation of a trajectory is proposed.


Introduction
The article deals with an optimal control problem with a fixed endpoint. The optimal control problem can be formulated as follows. Consider a controlled dynamical system where t ∈ [t 0 , T ], u = u(t) is vector control function belonging to some class of functions, for example, to the class of piecewise continuous functions taking values in some set: u(t) ∈ U ⊂ R k . Here and further, assume that the region of admissible states of the system Ω and function f (t, x, u) are such that solution of the Cauchy problem (1) uniquely exists in all cases under consideration.
Without loss of generality, consider the Mayer problem of the optimal control, as the Lagrange and the Bolza problems can be reduced to the Mayer problem by introducing an additional dynamic variable [Pontryagin et al., 1962]. Assume that the function f (t, x, u) is twice continuously differentiable over x and u. The Mayer problem consists in finding a control u(t) providing minimal value for a cost functional depending on x(T ) : The optimal control theory is closely related to the calculus of variations, from which it originates. For example, the principle of stationary action for a moving particle can be formulated as an optimal control problem, if the action is regarded as a cost functional, and the particle velocity is regarded as a control. In both of these areas, the first variation of trajectory of a dynamical system is commonly considered in many problems.
Particularly, various variational methods are developed for finding solution of an optimal control problem [Pytlak, 1999]. These methods use the first variation of the trajectory or the conjugate momentum Ψ(t) satisfying to the equation which is adjoint to the equation for the first variation, and to the final condition Orders of these methods are defined by orders of variation of the cost functional. If the first variation of the cost functional is taken into account only, such method is regarded as a first order methods. A method accounting also the second variation of the cost functional is a second order method. The commonly used expression for the second variation of the cost functional [Longmuir and Bohn, 1969;Gabasov and Kirillova, 1973;Miele, 1975;Golfetto and Silva, 2012] contains the first variation of the trajectory δx and the variation of the control function δu: Here index T denotes trabsposition, and H = Ψf. Here and further, we assume that the function g(x) in the expression (2) is twice continuously differentiable. When the control function is parametrized, this expression allows to find second derivations of the cost functional over the parameters used in process of computing of numerical solution. The second derivatives can be expressed through so the called 'matrix momenta' [Gabasov and Kirillova, 1973;Gorbunov and Lutoshkin, 2004], firstly introduced by G. Gabasov and F.M. Kirillova [Gabasov and Kirillova, 1973] for degenerate problems of optimal control. The matrix momenta represent components of some matrix satisfying to a corresponding boundary problems. Number of the boundary problems for the matrix momenta to be solved numerically is quadratic over number of the parameters. So far, in problems with a great number of parameters such approach is unrealizable. Therefore, it is reasonable try to use another expression for the second variation of the cost functional. In this article, we offer a new approach, where the second variation of the cost functional depends on the second variation of the trajectory. The main goals of this article is analysis of the second variation of the trajectory and application of it in numerical methods of finding solutions for the optimal control problem. The second variation of the trajectory is considered in the next section. Differential equation for the second variation is written. An integral representation for solution of this equation is obtained. The integral representation includes some tensor of the third rank and the Green function of differential equation for the first variation of the trajectory.
In the last section, a new method for numerical solution of the optimal control problem based on expression for the second variation of the cost functional including the second variation the trajectory is presented.

Second Variation of Trajectory
Let initial conditions or control function values are changed by small values of an order O(λ), where λ is small number. Assume that the disturbed solution can be represented in the form where x is solution of the undisturbed equation (1), δx is linear in λ, and δ 2 x is quadratic in λ. Let us call δ 2 x the second variation of the trajectory. We have Let us explain meaning of the notation used here and further. Consider the terms containing second derivatives of f. These terms represent tensors of the third rank. As polylinear forms, they have two arguments which are vectors, and one argument which is covector belonging to the conjugate space. Commonly, arguments from the conjugate space are being written on the left, and vector arguments consecutive on the right without separation by the comma. We will continue this tradition. For shortness, let denote two equal arguments by a square of one argument. For example, δx 2 should be regarded as δxδx. Besides, to emphasize that some expression is a vector argument, we will place it in parentheses if it consists from several cofactors. Such form of notation allows to reduce bulkiness of the following expressions.
Further, terms in right hand side of (7) contain tensor convolution over corresponding upper and lower indices: Upper and low indices are regarded as contravariant and covariant indices correspondingly. Further we will follow to the Einstein summation convention according to which summation is meant over all allowed values of repeated upper and lower indices. Substituting (6) and (7) into the equation (1) and equating terms of the first and of the second order in both side of this equation separately, we obtain two equations for terms of these orders: Here Solution of the linear nonuniform equation (8) can be written in the form where G(t, t ′ ) is the Green function satisfying to the equations and to the condition G(t, t) = E. Here, E is the tensor which components forms the identity matrix: and δ i j are Kronecker symbols. The equation (9) is also linear uniform equation for δ 2 x. Therefore, its solution can be written in the form Consider first term of the expression (14). It describes effect of interaction of two linear deviation of the solution due to nonlinearity of the equation (1). It results in appearing of a term of the second order. Write it in the form Integration in (15) is fulfilled over the domain Changing the order of integrations, write the expression (15) in the form Denote the integrand in the expression (16) by F (t ′ , t ′′ ) : It is obviously that Then the expression (16) can be written in the form Taking into account that Let introduce the tensor Then the expression (15) takes the form The tensor D is the tensor of the third rank, one time contravariant, and twice covariant, as the tensor ∂ 2 f /∂x 2 . This means that it has two vector arguments as a polylinear form. In the expression (22), they are written consequently on right, as it was pointed earlier. Note that the tensor D is symmetric relative to this vector arguments. It follows from its definition that D(t, t) = 0. Differentiating the expression (21) over the second argument being denoted here by t ′ , we get the following differential equation Here ellipsis denotes that one from the vector arguments of D is absent: the second one and the first one in the first term and in the second term of the right hand side of the (23) correspondingly. For greater clarity, let us write components of the expression for the second variation (14), of the tensor D, and of the tensor equation (23): 3 Second Order Method of Numerical Solution of the Optimal Control Problem Write change of the functional (2) under variation of the control u in the form Keeping in this equality linear terms, we have Introducing conjugate momentum satisfying to the adjoint equation (3) and the terminal condition (4), and taking into account the expression (13), we get for the first variation of the cost functional δΦ the well known expression Keeping in the equality (27) quadratic terms, we obtain the expression for second variation of the functional under variation of the control function The equality (31) can be written component-wise in the following form Let parameterize the control as follows. Let partition the segment [t 0 , T ] into subsegments [t i , t i+1 ], i = 1, M , t M = T, and let approximate the control by a piecewise constant vector function: The first derivatives of the functional over the parameters are determining on the base of the expression (30) as follows Analogically, making use the expression (32), one can find second derivatives of the functional over the parameters. Let us represent this expression in the form of sum of two terms Here the first term and the second term are contributions to the second derivative from the first term and from the second term of the expression (32) correspondingly. Let i > k. Then If i < k, than analogous expression can be obtained by permutation of indices in the expression for the previous case. Finally, if i = k, than For the second term, we have The algorithm of numerical solutions of the optimal control problem in the framework of presented approach consists in following steps. Firstly, some initial control piecewise vector function is chosen. At the next steps, it is being improved till some condition of termination of the calculations is met. At the next steps, control parameters are being changed. The vector of step in space of the parameters is being chosen in accordance with some method of optimization of the second order, that is a method which takes into account the first and the second derivatives of a function to be minimized.

Conclusion
The second variation of the trajectory of a dynamical system was considered also in the work [Henrion, 1975], where the equation for the second variation analogous to the equation (9) was written. In the present work, we found solution of this differential equation in integral form. The expression for the solution contains the Green function of equation for the first variation and the tensor D introduced here. The tensor D satisfies to the differential equation (23) and to some initial condition. For numerical solution of the optimal control problem, one should integrate all components of this tensor. The number of its components is cube of dimension of the phase space that is sufficiently smaller than in methods using the matrix momenta. By this reason we expect that the numerical methods based on using of the second variation of the trajectory will be effective in problems that require a great number of calculations, and when the control function can be parametrized by a great number of parameters. The optimal control theory has numerous physical applications. The second order method presented here can be applied for numerical solution of various physical and engineering problems. For example, it allows to find optimal trajectories for nonlinear dynamical systems, particularly, for aircraft [Miele, 1958;Gorbunov and Lutoshkin, 2004], spacecraft [Miele, 1958;Breakwell, 1962;Golfetto and Silva, 2012], and robots. Among other problems where this method can be used, it should be mentioned the problem of optimization of charged particle beam accelerator channel [Ovsyannikov, 1980;Bublik, Garashchenko, and Kirichenko, 1985;Drivotin et al., 1998;Ovsyannikov and Drivotin, 2003;Drivotin and Vlasova, 2014;Altsybeev et al.;Drivotin, 2018]. This problem was formulated as the problem of control of a dynamical system ensemble in the work [Ovsyannikov, 1980], and the numerical method based on the first variation of trajectory was developed [Ovsyannikov, 1980]. The method using the first variation requires a great amount of computation, espesially in this problem. Such method allows to compute gradient of the cost functional over parameters of the control function. From a geometric point of view, it is difficult to find direction of the descent in multidimensional space of parameters using only the gradient, because the gradient defines only a tangent plane to the surface on which functional value remains constant. To reduce amount of computations, it was proposed to apply for numerical solution the second order method based on the matrix momenta [Bublik, Garashchenko, and Kirichenko, 1985]. As it was explained previously, the method presented in this article is expected to be more effective than method using the matrix momenta. Although, this article is focused on application in numerical methods, we attach the highest importance to the second variation itself and hope that the integral representation for it will be useful in many problems of the optimal control.