Adaptive Synthesized Control for Solving The Optimal Control Problem

1. Introduction

Long ago Leonard Euler spoke about the optimal arrangement of everything in the world: "For since the fabric of the universe is most perfect and the work of a most wise Creator, nothing at all takes place in the universe in which some rule of maximum or minimum does not appear." Striving for optimality is natural in every sphere.

In order to optimally move an autonomous robot to a certain target position, currently, as a standard, engineers first solve the problem of optimal control, obtain the optimal trajectory, and then solve the additional problem of moving the robot along the obtained optimal trajectory. In most cases, the following approach is used to move the robot along a path. Initially the object is made stable relative to a certain point in the state space. Then the stability points are positioned along the desired path and the object is moved along the trajectory by switching these points from one point to another [1,2,3,4,5,6,7]. The difference between the existing methods is in solving the control synthesis problem to ensure stability relatively to some equilibrium point in the state space and in the location of these stability points.

Often, to ensure stability, the model of the control object is linearized relative to a certain point in the state space. Then, for the linear model of the object, a linear feedback control is found to arrange the eigenvalues of the closed-loop control system matrix on the left side of the complex plane. Sometimes, to improve the quality of stabilization, control channels or components of the control vector are defined that affect the movement of an object along a specific coordinate system axis of the state space. Then controllers, as a rule, PI controllers, are inserted into these channels with the coefficients that are adjusted according to the specified control quality criterion [3,4]. In some cases analytical or semi analytical methods are used to solve the control synthesis problem and building nonlinear stable control systems [5,7]. But the stability property of the nonlinear model of the control object, obtained from the linearization of this model, is generally preserved only in the vicinity of a stable equilibrium point.

The main drawback of this approach when the control object is moved along the stable points on the trajectory is that even if this trajectory is obtained as a result of solving the optimal control problem [8], then the movement itself will never be optimal. To ensure optimality, it is necessary to move along the trajectory at a certain speed, but when approaching the stable equilibrium point, the speed of the control object tends to zero.

The optimal control problem generally does not require ensuring the stability of the control object. The construction of a stabilization system that ensures the stability of the object relative to the equilibrium point in the state space is carried out by the researcher to achieve predictable behavior of the control object in the vicinity of a given trajectory.

The optimal control problem in the classical formulation is solved for a control object without any stabilization system, therefore, the resulting optimal control and the optimal trajectory will not be optimal for this object with further introduced stabilization system. It follows that the classical formulation of the optimal control problem [9] is missing something as far as its solution cannot be directly implemented in the real object, since this leads to an open-loop control. The open-loop control system is very sensitive to small disturbances, but they are always possible in real conditions, since no model accurately describes the control object. In order to achieve the optimal control in real object, it is necessary to build a feedback control system, which should provide some additional properties, for example stability relative to the trajectory or points on this trajectory. Authors [10,11] proposed an extended formulation statement of the optimal control problem, that has additional requirements established for the optimal trajectory. The optimal trajectory must have a non-empty neighborhood with a property of attraction. Performing these requirements provides implementation of the solution of the optimal control problem directly in the real control object.

In the work [12,13] an approach to solving the extended optimal control problem on the base of the synthesized control is presented. This approach ensures obtaining a solution of the optimal control problem in the class of practically implemented control functions. According to this approach, initially the control synthesis problem is solved. So the control object becomes stable in the state space relatively to some equilibrium point. In the second stage the optimal control problem is solved by determination of optimal positions of the stable equilibrium point. Switching stable points after a constant time interval ensures moving the control object from initial state to the terminal one optimally according to the given quality criterion. Optimal positions of stable equilibrium points can be far from the optimal trajectory in the state space, therefore a control object doesn’t slow down motion speed. Studies of synthesized control in various optimal control problems have shown that such control is not sensitive to perturbations and can be directly implemented in a real object [14,15].

In synthesized control the optimal control problem is solved for a control object already with a stabilization system. Other advantage of the synthesized control is a position of the stable point not changes in time interval, that is an optimal control function is solved in the class of piecewise constant functions, that simplifies a search of the optimal solution.

It is possible that piecewise constant control in the synthesized approach finds several optimal solutions with practically the same values of the quality criterion. This circumstance prompted us the idea to find among all almost optimal solutions one that is less sensitive to perturbations. This approach is called adaptive synthesized control.

In this work a principle of adaptive synthesized control is proposed in section 2, methods for solving it are discussed in section 3 and further in the section 4 a computational experiment of the solution of the optimal control problem for spatial motion of quadcopter by adaptive synthesized control is considered.

2. Adaptive Synthesized Control

Consider the principle of adaptive synthesized control for solving the optimal control problem in its extended formulation [10].

Initially, the control synthesis problem is solved to provide stability of control object relatively some point in the state space.

In the problem the mathematical model of the control object in the form of ordinary differential equation system is given.

$$\dot{\mathbf{x}}=\mathbf{f}(\mathbf{x},\mathbf{u}),$$

(1)

where $\mathbf{x}$ is a state vector, $\mathbf{x}\in {\mathbb{R}}^n$, $\mathbf{u}$ is a control vector, $\mathbf{u}\in \mathrm{U}\subseteq {\mathbb{R}}^m$, $\mathrm{U}$ is a compact set, that determines restrictions on the control vector.

The domain of admissible initial states is given

$${\mathrm{X}}_0\subseteq {\mathbb{R}}^n.$$

(2)

To solve the problem numerically, the initial domain (2) is taken in the form of the finite number of points in the state space

$${\tilde{\mathrm{X}}}_0=\{{\mathbf{x}}^{0,1},\dots ,{\mathbf{x}}^{0,K}\}.$$

(3)

Sometimes it is convenient to set one initial state and deviations from it:

$${\mathbf{x}}^{0,j}={\mathbf{x}}^0-{\Delta }_0+2\odot {\left(j\right)}_2{\Delta }_0,j=1,\dots ,2^n-1,$$

(4)

where ${\mathbf{x}}^0$ is a given initial state, ${\Delta }_0$ is a deviations vector, ${\Delta }_0={[{\Delta }_1\dots {\Delta }_n]}^T$, ⊙ is Hadamard product of vectors, ${\left(j\right)}_2$ is a binary code of the number j. In this case $K=2^n-1$.

The stabilization point as a terminal state is given

$${\mathbf{x}}^{f_1}\in {\mathbb{R}}^n.$$

(5)

It is necessary to find a control function in the form

$$\mathbf{u}=\mathbf{h}\left({\mathbf{x}}^{f_1}\mathbf{x}\right)\in \mathrm{U},$$

(6)

where $\mathbf{h}\left(\mathbf{x}\right):{\mathbb{R}}^n\to {\mathbb{R}}^m$, such that minimizes the quality criterion

$$J_0=\sum _{i=0}^K\left(t_{f_1,i}+p\parallel {\mathbf{x}}^{f_1}-\mathbf{x}(t_{f_1,i},{\mathbf{x}}^{0,i})\parallel \right)\to min,$$

(7)

where $t_{f_1,i}$ is a time of achieving the terminal state (5) from the initial state ${\mathbf{x}}^{0,i}$, $t_{f_1,i}$ is determined by an equation

$$t_{f_1,i}=\left\{\begin{array}{c}t,\mathrm{if}t<t^{+}\mathrm{and}\parallel {\mathbf{x}}^{f_1}-\mathbf{x}(t,{\mathbf{x}}^{0,i})\parallel \le {\epsilon }_0\hfill \\ t^{+},\mathrm{otherwise}\hfill \end{array}\right.,$$

(8)

$\mathbf{x}(t,{\mathbf{x}}^{0,i})$ is a particular solution from initial state ${\mathbf{x}}^{0,i}$, $i=1,\dots ,K$, of the differential equation (1) with inserted there control function (6)

$$\dot{\mathbf{x}}=\mathbf{f}(\mathbf{x},\mathbf{h}({\mathbf{x}}^{f_1}-\mathbf{x})),$$

(9)

${\epsilon }_0$ is a given accuracy for hitting to terminal state (5), $t^{+}$ is a given maximal time for control process, p is a weight coefficient.

Further in the principle of the synthesized optimal control the following optimal control problem is considered. The model of control object in the form (9) is used

$$\dot{\mathbf{x}}=\mathbf{f}(\mathbf{x},\mathbf{h}({\mathbf{x}}^{*}-\mathbf{x})),$$

(10)

where the terminal state vector (5) is changed into the new unknown vector ${\mathbf{x}}^{*}$, that will be a control vector in the considered optimal control problem.

In accordance with the classical formulation of the optimal control problem, the initial state of the object (10) is given

$${\mathbf{x}}^0\in {\mathbb{R}}^n.$$

(11)

As far as in the engineering practice there can be some deviations in the initial position, therefore, in the adaptive synthesized control instead of one initial state (11) the set of initial states are used defined by the equation (4). The vector of initial deviations ${\Delta }_0$ is defined as a level of disturbances.

The goal of control is defined by achievement of the terminal state

$${\mathbf{x}}^f\in {\mathbb{R}}^n.$$

(12)

The quality criterion is given

$$J_1=\underset{0}{\overset{t_f}{\int }}{f}_0(\mathbf{x},\mathbf{h}({\mathbf{x}}^{*}-\mathbf{x}))dt\to min,$$

(13)

where $t_f$ is a terminal time, $t_f$ is not given but is limited, $t_f\le t^{+}$, $t^{+}$ is a given limit time of control process.

According to the principle of synthesized control, it is necessary to choose time interval $\Delta t$ and to search for optimal constant values of the control vector ${\mathbf{x}}^{*,i}$ for each interval

$${\mathbf{x}}^{*}={\mathbf{x}}^{*,i},\mathrm{if}(i-1)\Delta t\le t<i\Delta t,i=1\dots ,M,$$

(14)

where M is a number of intervals

$$M=⌊\frac{t^{+}}{\Delta t}⌋,$$

(15)

such, that particular solution of the system (10) from the given initial state (11) reaches terminal state with optimal value of the quality criterion (13).

Algorithmically in the second stage of the adaptive synthesized control approach, the optimal values of the vector ${\mathbf{x}}^{*}$ are found as a result of the optimization task with the following quality criterion, which takes into account the given grid according to the initial conditions

$$J_2=\sum _{i=1}^K\left(\underset{0}{\overset{t_{f,i}}{\int }}{f}_0(\mathbf{x},\mathbf{h}({\mathbf{x}}^{*,i}-\mathbf{x}))dt+p\parallel {\mathbf{x}}^f-\mathbf{x}(t_{f,i},{\mathbf{x}}^{0,i})\parallel \right)\to \underset{{\mathbf{x}}^{*}}{min},$$

(16)

where K is number of initial states, $t_{f,i}$ is determined by equation (8).

3. Methods of Solving

As described in the previous section, the approach based on the principle of adaptive synthesized optimal control consists of two stages.

To implement the first stage of the approach under consideration for solving the control synthesis problem (1)–(9) any known method can be used mostly depending on the type of the model used to describe the control object. In practice, when developing a stabilization system, engineers do not use the initial performance criterion (7), since the main goal is to ensure stability. To solve the synthesis problem and obtain an equilibrium point, methods of modal control [16] can be applied for linear systems, as well as other analytical methods such as backstepping [17] or synthesis based on the application of the Lyapunov function [18]. In practice, stability is ensured through linearization of the model (1) in the terminal state and setting PI or PID-controllers in control channels [19,20]. The mathematical formulation of the stabilization problem as a control synthesis problem is needed to apply numerical methods and automatically obtain a feedback control function. Today, to solve the synthesis problem for nonlinear dynamic objects of varying complexity, modern numerical methods of machine learning can be applied [21]. In the present paper machine learning by symbolic regression [22] is used.

Symbolic regression allows to find a mathematical expression for desired function in the form of special code. All symbolic regression methods differs in the form of this code. The search for solutions is performed in the space of codes by a special genetic algorithm.

Let us demonstrate main features of symbolic regression on the example of the network operator method (NOP), that was used in this work in the computational experiment. To code a mathematical expression NOP uses an alphabet of elementary functions:

–

functions without arguments or parameters and variables of the mathematical expression;

$${\mathrm{F}}_0=\{f_{0,1}=x_1,\dots ,f_{0,n}=x_n,f_{0,n+1}=q_1,\dots ,f_{0,n+r}=q_{n+r}\};$$

(17)

–

functions with one argument

$${\mathrm{F}}_1=\{f_{1,1}\left(z\right)=z,f_{1,2}\left(z\right),\dots ,f_{1,W}\left(z\right)\};$$

(18)

–

function with two arguments

$${\mathrm{F}}_2=\{f_{2,1}(z_1,z_2),\dots ,f_{2,V}(z_1,z_2)\}.$$

(19)

Any elementary function is coded by two digits, the first one is the number of arguments, the second one is the function number in corresponding set. These digits are written as indexes of elements in the introduced sets of the alphabet (17)–(19). The set of functions with one argument must include the identity function $f_{1,1}\left(z\right)=1$. Functions with two arguments should be commutative, associative and have a unit element.

NOP encodes a mathematical expression in the form of an oriented graph. Source-nodes of the NOP-graph are connected with functions without arguments, other nodes are connected with functions with two arguments. Arcs of the NOP-graph are connected with functions with one argument. If on the NOP-graph some node has one input arc, then the second argument is a unit element for the function with two arguments connected with this node.

Let us define the following alphabet of elementary functions:

$$\begin{array}{ccc}\hfill {\mathrm{F}}_0& =& \{f_{0,1}=x_1,f_{0,2}=x_2,f_{0,3}=q_1,f_{0,4}=q_2\},\hfill \\ \hfill {\mathrm{F}}_1& =& \{f_{1,1}\left(z\right)=z,f_{1,2}\left(z\right)=-z,f_{1,3}\left(z\right)=cos\left(z\right),f_{1,4}\left(z\right)=sin\left(z\right)\},\hfill \\ \hfill {\mathrm{F}}_2& =& \{f_{2,1}(z_1,z_2)=z_1+z_2,f_{2,2}(z_1,z_2)=z_1{z}_2\}.\hfill \end{array}$$

(20)

With this alphabet the following mathematical expressions can be encoded in the form of NOP:

$$\begin{array}{ccc}\hfill y_1& =& cos\left(x_1\right)sin\left(x_2\right)-sin\left(x_1\right)cos\left(x_2\right),\hfill \\ \hfill y_2& =& cos(q_1{x}_1+q_{2}sin\left(x_2\right)),\hfill \\ \hfill y_3& =& q_{1}sin\left(q_2{x}_2\right)+q_{2}cos\left(q_1{x}_1\right),\hfill \\ \hfill y_4& =& q_{1}sin(q_{2}cos\left(x_1\right))+q_{2}cos(q_{1}sin\left(x_2\right)).\hfill \end{array}$$

(21)

The NOP-graphs of these mathematical expressions are presented in Figure 1. In the computer memory the NOP-graphs are presented in the form of integer matrices.

$$\mathbf{\Psi }=\left[{\psi }_{i,j}\right],i,j=1,\dots ,L.$$

(22)

As far as the NOP-nodes are enumerated in such a way that the node number from which an arc comes out is less than the node number to which an arc enters, then the NOP matrix has an upper triangular form. Every line of the matrix corresponds some node of the graph. Lines with zeros in the main diagonal corresponds to source-nodes of the graph. Other elements in the main diagonal are the function numbers with two arguments. Non zero elements above the main diagonal are the function numbers with one argument.

NOP matrices for the mathematical expressions (21) have the following forms

$${\mathbf{\Psi }}_1=\left[\begin{array}{ccccc}0& 0& 3& 4& 0\\ 0& 0& 4& 3& 0\\ 0& 0& 2& 0& 2\\ 0& 0& 0& 2& 1\\ 0& 0& 0& 0& 1\end{array}\right],{\mathbf{\Psi }}_2=\left[\begin{array}{cccccccc}0& 0& 0& 0& 1& 0& 0& 0\\ 0& 0& 0& 0& 0& 4& 0& 0\\ 0& 0& 0& 0& 1& 0& 0& 0\\ 0& 0& 0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& 2& 0& 1& 0\\ 0& 0& 0& 0& 0& 2& 1& 0\\ 0& 0& 0& 0& 0& 0& 1& 3\\ 0& 0& 0& 0& 0& 0& 0& 1\end{array}\right],$$

$${\mathbf{\Psi }}_3=\left[\begin{array}{ccccccccc}0& 0& 0& 0& 1& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 1& 0& 0& 0\\ 0& 0& 0& 0& 1& 0& 0& 1& 0\\ 0& 0& 0& 0& 0& 1& 1& 0& 0\\ 0& 0& 0& 0& 2& 0& 3& 0& 0\\ 0& 0& 0& 0& 0& 2& 0& 4& 0\\ 0& 0& 0& 0& 0& 0& 2& 0& 1\\ 0& 0& 0& 0& 0& 0& 0& 2& 1\\ 0& 0& 0& 0& 0& 0& 0& 0& 1\end{array}\right],{\mathbf{\Psi }}_4=\left[\begin{array}{ccccccccc}0& 0& 0& 0& 3& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 4& 0& 0& 0\\ 0& 0& 0& 0& 0& 1& 1& 0& 0\\ 0& 0& 0& 0& 1& 0& 0& 1& 0\\ 0& 0& 0& 0& 2& 0& 4& 0& 0\\ 0& 0& 0& 0& 0& 2& 0& 3& 0\\ 0& 0& 0& 0& 0& 0& 2& 0& 1\\ 0& 0& 0& 0& 0& 0& 0& 2& 1\\ 0& 0& 0& 0& 0& 0& 0& 0& 1\end{array}\right].$$

(23)

To calculate a mathematical expression by its NOP matrix, initially, the vector of nodes is determined. Number of components of the vector of nodes equals to the number of nodes in a graph. The initial vector of nodes includes variables and parameters in positions that correspond to source nodes, as well as other components equal to the unit elements of the corresponding functions with two arguments. Further every line of matrix is checked. If element of matrix doesn’t equal to zero, then corresponding element of the vector of nodes is changed. To calculate mathematical expression by the NOP matrix the following equation is used

$$z_j^{\left(i\right)}\leftarrow \left\{\begin{array}{c}{f}_{2,{\psi }_{j,j}}(z_j^{(i-1)},f_{1,{\psi }_{i,j}}\left(z^{{(i-1)}_i}\right)),\mathrm{if}{\psi }_{i,j}\ne 0\hfill \\ z_j^{(i-1)},\mathrm{otherwise}\hfill \end{array}\right.,i=1,\dots ,L-1,j=i+1,\dots ,L,$$

(24)

where

$$z_i^{\left(0\right)}=\left\{\begin{array}{c}{f}_{0,i},\mathrm{if}{\psi }_{i,i}=0\hfill \\ e_{{\psi }_{i,i}},\mathrm{otherwise}\hfill \end{array}\right.,$$

(25)

$e_j$ is a unit element for function with two arguments $f_{2,j}(z_1,z_2)$,

$$f_{2,j}(e_j,z)=f(z,e_j)=z.$$

(26)

Consider an example of calculating the second mathematical expression in (21) on its NOP matrix ${\mathbf{\Psi }}_2$.

Initial vector of nodes is

$${\mathbf{z}}^{\left(0\right)}={\left[x_1{x}_2{q}_1{q}_{2}1100\right]}^T.$$

Further all strings in the matrix ${\mathbf{\Psi }}_2$ are checked and non zero elements are found.

$${\psi }_{1,5}=1,z_5^{\left(1\right)}=f_{2,2}(z_5^{\left(0\right)},f_{1,1}\left(z_1^{\left(0\right)}\right))=1·f_{1,1}\left(z_1^{\left(0\right)}\right)=1·x_1=x_1,$$

$${\psi }_{2,6}=4,z_6^{\left(2\right)}=f_{2,2}(z_6^{\left(1\right)},f_{1,4}\left(z_1^{\left(1\right)}\right))=1·f_{1,4}\left(z^{\left(1\right)}\right)=1·sin\left(x_2\right)=sin\left(x_2\right),$$

$${\psi }_{3,5}=1,z_5^{\left(3\right)}=f_{2,2}(z_5^{\left(2\right)},f_{1,1}\left(z_3^{\left(2\right)}\right))=x_1·q_1=q_1{x}_1,$$

$${\psi }_{4,6}=1,z_6^{\left(4\right)}=f_{2,2}(z_6^{\left(3\right)},f_{1,1}\left(z_4^{\left(3\right)}\right))=sin\left(x_2\right)·q_2=q_{2}sin\left(x_2\right),$$

$${\psi }_{5,7}=1,z_7^{\left(5\right)}=f_{2,1}(z_7^{\left(4\right)},f_{1,1}\left(z_5^{\left(4\right)}\right))=0+q_1{x}_1=q_1{x}_1,$$

$${\psi }_{6,7}=1,z_7^{\left(6\right)}=f_{2,1}(z_7^{\left(5\right)},f_{1,1}\left(z_6^{\left(5\right)}\right))=q_1{x}_1+q_{2}sin\left(x_2\right),$$

$${\psi }_{7,8}=3,z_8^{\left(7\right)}=f_{2,1}(z_8^{\left(6\right)},f_{1,3}\left(z_7^{\left(6\right)}\right))=0+cos(q_1{x}_1+q_{2}sin\left(x_2\right))=cos(q_1{x}_1+q_{2}sin\left(x_2\right)).$$

The last mathematical expression coincides with the needed mathematical expression for $y_2$ (21).

So, we considered the way of codding in the NOP method. Then, to search for an optimal mathematical expression in some task, the NOP method applies a principle of small variations of a basic solution. According to this principle one possible solution is encoded in the form of the NOP matrix ${\mathbf{\Psi }}_0$. This solution is the basic solution and it is set by a researcher as some good solution. Others possible solutions are presented in the form of sets of small variation vectors. A small variation vector consists of four integer numbers

$$\mathbf{w}={\left[w_1{w}_2{w}_3{w}_4\right]}^T,$$

(27)

where $w_1$ is a type of small variation, $w_2$ the a line number of NOP-matrix, $w_3$ is a column number of NOP-matrix, $w_4$ is a new value of a NOP-matrix element. There are four types of small variations: $w_1=0$ is an exchange of the function with one argument: if ${\psi }_{w_2,w_3}\ne 0$, then ${\psi }_{w_2,w_3}\leftarrow w_4$; $w_1=1$ is an exchange of the function with two arguments: if ${\psi }_{w_2,w_2}\ne 0$, then ${\psi }_{w_2,w_2}\leftarrow w_4$; $w_1=2$ is an insertion of the additional function with one argument: if ${\psi }_{w_2,w_3}=0$, then ${\psi }_{w_2,w_3}\leftarrow w_4$; $w_1=3$ is an elimination of the function with one argument: if ${\psi }_{w_2,w_3}\ne 0$ and $\exists {\psi }_{w_2,j}\ne 0$, $j>w_2$, $j\ne w_3$ and $\exists {\psi }_{i,w_3}\ne 0$, $i\ne w_2$, then ${\psi }_{w_2,w_3}\leftarrow 0$.

The initial population includes H possible solutions. Each possible solution $i\in \{1,\dots ,H\}$ except the basic solution is encoded in the form of the set of small variation vectors

$${\mathrm{W}}_i=({\mathbf{w}}^{i,1,},\dots ,{\mathbf{w}}^{i,d}),i\in \{1,\dots ,H\},$$

(28)

where d is a depth of variations, which is set as a parameter of the algorithm.

The NOP-matrix of a possible solution is determined after application of all small variations to the basic solution

$${\mathbf{\Psi }}_i={\mathbf{w}}^{i,d}\circ \dots \circ {\mathbf{w}}^{i,1}\circ {\mathbf{\Psi }}_0,i\in \{1,\dots ,H\},$$

(29)

here small variation vector is written as a mathematical operator changing matrix ${\mathbf{\Psi }}_0$.

During the search process sometimes the basic solution is replaced by the current best possible solution. This process is called a change of an epoch.

Consider an example of applying small variations to the NOP-matrix ${\mathbf{\Psi }}_3$. Let $d=3$ and there are three following small variation vectors

$${\mathbf{w}}^1={\left[0352\right]}^T,{\mathbf{w}}^2={\left[2563\right]}^T,{\mathbf{w}}^3={\left[2894\right]}^T$$

After application of these small variation vectors to the NOP-matrix ${\mathbf{\Psi }}_3$ the following NOP-matrix is obtained

$${\mathbf{\Psi }}_5=\left[\begin{array}{ccccccccc}0& 0& 0& 0& 1& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 1& 0& 0& 0\\ 0& 0& 0& 0& 2& 0& 0& 1& 0\\ 0& 0& 0& 0& 0& 1& 1& 0& 0\\ 0& 0& 0& 0& 2& 3& 3& 0& 0\\ 0& 0& 0& 0& 0& 2& 0& 4& 0\\ 0& 0& 0& 0& 0& 0& 2& 0& 1\\ 0& 0& 0& 0& 0& 0& 0& 2& 4\\ 0& 0& 0& 0& 0& 0& 0& 0& 1\end{array}\right].$$

This NOP-matrix corresponds to the following mathematical expression

$$y_5=sin(q_{1}sin(q_2{x}_{2}cos(-q_1{x}_1)))+q_{2}cos(-q_1{x}_1).$$

As a search engine for the optimal solution a genetic algorithm is used. To perform a crossover operation two possible solutions are selected randomly

$$\begin{array}{ccc}\hfill {\mathrm{W}}_{\alpha }& =& ({\mathbf{w}}^{\alpha ,1},\dots ,{\mathbf{w}}^{\alpha ,d}),\hfill \\ \hfill {\mathrm{W}}_{\beta }& =& ({\mathbf{w}}^{\beta ,1},\dots ,{\mathbf{w}}^{\beta ,d}).\hfill \end{array}$$

(30)

A crossover point is selected randomly $c\in \{1,\dots ,d\}$. Two new possible solutions are obtained as the result of exchanging elements of the selected possible solutions after crossover point:

$$\begin{array}{ccc}\hfill {\mathrm{W}}_{H+1}& =& ({\mathbf{w}}^{\alpha ,1},\dots ,{\mathbf{w}}^{\alpha ,c},{\mathbf{w}}^{\beta ,c+1},\dots ,{\mathbf{w}}^{\beta ,d}),\hfill \\ \hfill {\mathrm{W}}_{H+2}& =& ({\mathbf{w}}^{\beta ,1},\dots ,{\mathbf{w}}^{\beta ,c},{\mathbf{w}}^{\alpha ,c+1},\dots ,{\mathbf{w}}^{\alpha ,d}).\hfill \end{array}$$

(31)

The second stage of the synthesized principle under consideration is to solve the problem of optimal control via determination of the optimal position of the equilibrium points. Studies have shown [23] that for a complex optimal control problem with phase constraints, evolutionary algorithms allow to cope with such problems. Good results were demonstrated [24] by such algorithms as a genetic algorithm (GA) [25], a particle swarm optimization (PSO) algorithm [26], a grey wolf optimizer (GWO) algorithm [27] or a hybrid algorithm [28] involving one population of possible solutions and all three evolutionary transformations of GA, PSO and GWO selected randomly.

4. Computational Experiment

Consider the optimal control problem for spatial motion of a quadcopter. In the problem, the quadcopter should move in a minimum time on a closed loop circle from the given initial state to the the same terminal state, avoiding collisions with obstacles and passing through the given areas.

4.1. Mathematical Model of Spatial Movement of Quadcopter

In general case the mathematical model of quadcopter as a hard body has the following form:

$$\begin{array}{ccc}\hfill \ddot{x}& =& F(cos(\gamma )sin(\theta )cos(\psi )+sin(\gamma )sin(\psi \left)\right)/m,\hfill \\ \hfill \ddot{y}& =& Fcos\left(\gamma \right)cos\left(\theta \right)/m-g,\hfill \\ \hfill \ddot{z}& =& F(cos(\gamma )sin(\theta )sin(\psi )+sin(\gamma )cos(\psi \left)\right)/m,\hfill \\ \hfill \ddot{\gamma }& =& ((I_{yy}+I_{zz})\dot{\theta }\dot{\psi }+M_x)/I_{xx},\hfill \\ \hfill \ddot{\psi }& =& ((I_{zz}+I_{xx})\dot{\gamma }\dot{\theta }+M_y)/I_{yy},\hfill \\ \hfill \ddot{\theta }& =& ((I_{xx}+I_{yy})\dot{\gamma }\dot{\theta }+M_z/I_{zz},\hfill \end{array}$$

(32)

where F is a summary thrust force of all drone screws, m is a mass of drone, g is acceleration of gravity, $g=9.80665$, $M_x$, $M_y$, $M_z$ are control moments around of respective axes.

The Figure 2 shows how angles of quadcopter turn are linked with its axes.

To transform the model (32) to a vector record the following designations are entered: $x=x_1$, $y=x_2$, $z=x_3$, ${\dot{x}}_1=x_4$, ${\dot{x}}_2=x_5$, ${\dot{x}}_3=x_6$, $\gamma =x_7$, $\psi =x_8$, $\theta =x_9$, $\dot{\gamma }=x_{10}$, $\dot{\psi }=x_{11}$, $\dot{\theta }=x_{12}$, $M_1=M_x$, $M_2=M_y$, $M_3=M_z$.

As a result the following mathematical model is received:

$$\begin{array}{ccc}\hfill {\dot{x}}_1& =& x_4,\hfill \\ \hfill {\dot{x}}_2& =& x_5,\hfill \\ \hfill {\dot{x}}_3& =& x_6,\hfill \\ \hfill {\dot{x}}_4& =& F(cos\left(x_7\right)sin\left(x_9\right)cos\left(x_8\right)+sin\left(x_7\right)sin\left(x_8\right))/m,\hfill \\ \hfill {\dot{x}}_5& =& F(cos\left(x_7\right)cos\left(x_9\right)/m-g,\hfill \\ \hfill {\dot{x}}_6& =& F(cos\left(x_7\right)sin\left(x_9\right)sin\left(x_8\right)+sin\left(x_7\right)cos\left(x_8\right))/m,\hfill \\ \hfill {\dot{x}}_7& =& x_{10},\hfill \\ \hfill {\dot{x}}_8& =& x_{11},\hfill \\ \hfill {\dot{x}}_9& =& x_{12},\hfill \\ \hfill {\dot{x}}_{10}& =& ((I_{yy}+I_{zz})x_{11}{x}_{12}+M_x)/I_{xx},\hfill \\ \hfill {\dot{x}}_{11}& =& ((I_{zz}+I_{xx})x_{10}{x}_{12}+M_y)/I_{yy},\hfill \\ \hfill {\dot{x}}_{12}& =& ((I_{xx}+I_{yy})x_{10}{x}_{11}+M_z)/I_{zz},\hfill \end{array}$$

(33)

where $\mathbf{x}$ is a state space vector, $\mathbf{x}={[x_1\dots x_n]}^T$, $\mathbf{M}$ is a vector of control moments, $\mathbf{M}={\left[M_1{M}_2{M}_3\right]}^T$.

As a rule quadcopters are manifactured with some angle stabilization systems. It means that a drone can be stabilize on any angle from some interval. The system of angles stabilization provides stable location of drone relatively given angles by control moments

$$M_i=w_i(x_7^{*}-x_7,x_8^{*}-x_8,x_9^{*}-x_9,x_{10},x_{11},x_{12}),i=1,2,3.$$

(34)

Assume that the angular stabilization system works out the given angles of the quadcopter quickly enough, at least in comparison with spatial movement. In this case we can assume that the control of the spatial movement of the quadcopter is carried out using the angular position of the drone and the thrust force. Let’s define components of the spatial control vector: $x_7=u_1$, $x_8=u_2$, $x_9=u_3$, $F/m=u_4$.

In result we receive the following model of spatial quadcopter movement

$$\begin{array}{ccc}\hfill {\dot{x}}_1& =& x_4,\hfill \\ \hfill {\dot{x}}_2& =& x_5,\hfill \\ \hfill {\dot{x}}_3& =& x_6,\hfill \\ \hfill {\dot{x}}_4& =& u_4(cos\left(u_1\right)sin\left(u_3\right)cos\left(u_2\right)+sin\left(u_1\right)sin\left(u_2\right)),\hfill \\ \hfill {\dot{x}}_5& =& u_{4}cos\left(u_1\right)cos\left(u_3\right)-g,\hfill \\ \hfill {\dot{x}}_6& =& u_4(cos\left(u_1\right)sin\left(u_3\right)sin\left(u_2\right)+sin\left(u_1\right)cos\left(u_2\right)).\hfill \end{array}$$

(35)

In the work this model is used in the control problem of quadcopter spatial motion.

4.2. The optimal control problem for spatial motion of quadcopter

The model (35) of control object is given. Here $\mathbf{x}$ is a state space vector, $\mathbf{x}\in {\mathbb{R}}^6$, $\mathbf{u}$ is a control vector $\in \mathrm{U}\in {\mathbb{R}}^4$. $\mathrm{U}$ is a compact set, that defines restrictions on values of control vector components,

$$\begin{array}{ccccc}\hfill u_1^{-}=-\pi /12& \le & u_1& \le & \pi /12=u_1^{+},\hfill \\ \hfill u_2^{-}=-\pi & \le & u_2& \le & \pi =u_2^{+},\hfill \\ \hfill u_3^{-}=-\pi /12& \le & u_3& \le & \pi /12=u_3^{+},\hfill \\ \hfill u_4^{-}=0& \le & u_4& \le & 12=u_4^{+}.\hfill \end{array}$$

(36)

According to the principle of synthesized control, initially the control synthesized problem (1)–(9) is solved. The model (35) is used as a model of control object. To construct the set of initial states (4) the following vector of deviations is used

$${\Delta }_0={\left[222000\right]}^T.$$

(37)

In the problem initial state and terminal state were equal

$${\mathbf{x}}^0={\mathbf{x}}^f={\left[050000\right]}^T.$$

(38)

For calculation of the quality criterion (7) the following parameters are used $t^{+}=2$, ${\epsilon }_0=0.1$, $p=2$.

To solve the control synthesis problem the network operator method[22] was used. NOP found the following solution

$$u_i=\left\{\begin{array}{c}{u}_i^{+},\mathrm{if}{\widehat{u}}_i>u_i^{+}\hfill \\ u_i^{-},\mathrm{if}{\widehat{u}}_i<u_i^{-}\hfill \\ {\widehat{u}}_i,\mathrm{otherwise}\hfill \end{array}\right.,i=1,\dots ,m=4,$$

(39)

where

$${\widehat{u}}_1=\mu \left(C\right),$$

(40)

$${\widehat{u}}_2={\widehat{u}}_1-{\widehat{u}}_1^3,$$

(41)

$${\widehat{u}}_3={\widehat{u}}_2+{\rho }_{19}(W+\mu \left(C\right))+{\rho }_{17}\left(A\right),$$

(42)

$${\widehat{u}}_4={\widehat{u}}_3+ln\left(\right|{\widehat{u}}_2\left|\right)+\mathrm{sgn}(W+\mu \left(C\right))\sqrt{|W+\mu (C\left)\right|}+{\rho }_{19}\left(W\right)+$$

$$arctan\left(H\right)+\mathrm{sgn}\left(F\right)+arctan\left(E\right)+exp\left(q_2(x_2^f-x_2)\right)+\sqrt{q_1},$$

(43)

$$C=q_6(x_6^f-x_6)+q_3(x_3^f-x_3),W=V+tanh\left(G\right)+exp\left(D\right),$$

$$A=q_1(x_1^f-x_1)+q_4(x_4^f-x_4),H=G+tanh\left(F\right)+{\rho }_{18}\left(B\right),$$

$$F=E+C+arctan\left(D\right)-B,E=D+\mathrm{sgn}(x_5^f-x_5)+{(x_2^f-x_2)}^3,$$

$$V=exp\left(H\right)+cos\left(q_6(x_6^f-x_6)\right)+\mathrm{sgn}\left(D\right)\sqrt{\left|D\right|},G=F+\sqrt[3]{E}+sin\left(A\right),$$

$$B=sin\left(q_6(x_6^f-x_6)\right)+q_5(x_5^f-x_5)+q_2(x_2^f-x_2)+cos\left(q_1\right)+\vartheta (x_2^f-x_2),$$

$$D={\rho }_{17}\left(C\right)+B^3+A+\vartheta \left(q_5(x_5^f-x_5)\right)+{(x_5^f-x_5)}^2,$$

$$\mu \left(z\right)=\left\{\begin{array}{c}z,\mathrm{if}\left|z\right|<1\hfill \\ \mathrm{sgn}\left(z\right),\mathrm{otherwise}\hfill \end{array}\right.,{\rho }_{17}\left(z\right)=\mathrm{sgn}\left(z\right)ln\left(\right|z|+1),$$

$${\rho }_{18}\left(z\right)=\mathrm{sgn}\left(z\right)(exp(\left|z\right|)-1),{\rho }_{19}\left(z\right)=\mathrm{sgn}\left(z\right)exp(-|z\left|\right),$$

$q_1=7.26709$, $q_2=11.46143$, $q_3=12.77026$, $q_4=3.20630$, $q_5=8.38501$, $q_6=5.56250$.

In the second stage the optimal control problem is considered. In the problem the mathematical model (35) is given. The initial state coincides with terminal state (38).

It is necessary to find a control in the form of points in the states space (14). For synthesized control it is necessary to minimize the following quality criterion

$$J_3=t_f+p_1\parallel {\mathbf{x}}^f-\mathbf{x}+\left(t_f\right)\parallel +p_2\sum _{i=0}^N\underset{0}{\overset{t_f}{\int }}\vartheta \left({\phi }_i\left(\mathbf{x}\right)\right)dt+$$

$$p_3\sum _{j=1}^S{p}_3\vartheta (\underset{t}{min}|{\delta }_j\left(\mathbf{x}\right)|-\epsilon )\to \underset{{\mathbf{x}}^{*}}{min},$$

(44)

where $p_1=2$, $p_2=3$, $p_3=3$,

$${\phi }_i\left(\mathbf{x}\right)=r_i-\sqrt{{(x_{i,1}-x_1)}^2+{(x_{i,3}-x_3)}^2},$$

(45)

$i=1,\dots ,N=4$, $r_1=r_2=r_3=r_4=2$, $x_{1,1}=5$, $x_{1,3}=0$, $x_{2,1}=10$, $x_{2,3}=5$, $x_{3,1}=5$, $x_{3,3}=10$, $x_{4,1}=0$, $x_{4,3}=5$,

$${\delta }_j\left(\mathbf{x}\right)=\sqrt{{(y_{i,1}-x_1)}^2+{(y_{j,3}-x_3)}^2},$$

(46)

$j=1,\dots ,S=7$, $y_{1,1}=5$,- $y_{1,3}=-2$, $y_{2,1}=10$, $y_{2,3}=0$, $y_{3,1}=12$,+ $y_{3,3}=5$, $y_{4,1}=10$, $y_{4,3}=10$, $y_{5,1}=5$, $y_{5,3}=12$, $y_{6,1}=0$, $y_{6,3}=10$, $y_{7,1}=2$, $y_{7,3}=5$, $\epsilon =0.6$.

In the optimal control problem the terminal time $t_f$ is determined by the equation (8) with $t^{+}=14.4$, ${\epsilon }_0=0.1$. It is necessary to find coordinates of control points on each time interval, $\Delta t=0.8$. The desired vector includes $3M$ parameters, where

$$M=⌊\frac{t^{+}}{\Delta t}⌋=⌊\frac{14.4}{0.8}⌋=18,$$

(47)

that is ${\mathbf{q}}^{*}={[q_1\dots q_{54}]}^T$. The hybrid evolutionary algorithm has found the following optimal solution

$$\begin{array}{c}{\mathbf{x}}^{*,1}={[4.839101.14025-5.22899]}^T,{\mathbf{x}}^{*,2}={[11.070566.79389-2.48647]}^T,\\ {\mathbf{x}}^{*,3}={\left[9.198081.5467415.87195\right]}^T,{\mathbf{x}}^{*,4}={[-0.122040.12276-1.82381]}^T,\\ {\mathbf{x}}^{*,5}={[-4.083472.936585.89553]}^T,{\mathbf{x}}^{*,6}={\left[16.728962.180222.27907\right]}^T,\\ {\mathbf{x}}^{*,7}={\left[1.181062.5658214.41088\right]}^T,{\mathbf{x}}^{*,8}={[8.671985.78737-2.90409]}^T,\\ {\mathbf{x}}^{*,9}={\left[8.594782.7394811.33252\right]}^T,{\mathbf{x}}^{*,10}{=[-1.25924[-1.97448-1.42747]}^T,\\ {\mathbf{x}}^{*,11}={[2.454457.42257-0.38164]}^T,{\mathbf{x}}^{*,12}={[8.68306-0.7849615.41667]}^T,\\ {\mathbf{x}}^{*,13}={\left[0.609727.027247.66403\right]}^T,{\mathbf{x}}^{*,141}={[-0.599750.39324-1.31307]}^T,\\ {\mathbf{x}}^{*,15}={[-2.390047.952793.02003]}^T,{\mathbf{x}}^{*,16}={\left[2.526426.693329.17356\right]}^T\\ {\mathbf{x}}^{*,17}={[-0.958964.42529-0.36318]}^T,{\mathbf{x}}^{*,18}={[-0.011935.0282115.40007]}^T.\end{array}$$

(48)

For the found solution (48) the value of the quality criterion is $J_3=14.7010$.

In the Figure 3 projections of the optimal trajectory on the horizontal plane $\{x_1;x_3\}$ are presented. Here red circles are phase constraints described by the equation (45), small black circle are passing areas described by the equation (delta), small black squares are control points (48).

For the new adaptive synthesized control, proposed in this paper, the set of initial states is determined by the equation (3) with deviation vector

$${\Delta }_0={\left[0.20.20.2000\right]}^T.$$

(49)

It is necessary to find the same number of control points according to the following quality criterion

$$J_4=\sum _{k=1}^K\left(t_{f,k}+p_1\parallel {\mathbf{x}}^f-\mathbf{x}(t_{f,i},{\mathbf{x}}^{0,i})+p_2\sum _{i=0}^N\underset{0}{\overset{t_f}{\int }}\vartheta \left({\phi }_i\left(\mathbf{x}\right)\right)dt+\right.$$

$$\left.p_3\sum _{j=1}^S{p}_3\vartheta (\underset{t}{min}|{\delta }_j\left(\mathbf{x}\right)|-\epsilon )\right)\to \underset{{\mathbf{x}}^{*}}{min},$$

(50)

where $K=7$, $t_{f,k}$ is defined by equation (8). Other parameters of the criterion are the same as for the criterion (44).

The hybrid evolutionary algorithm has found the following optimal solution

$$\begin{array}{c}{\mathbf{x}}^{*,1}={[17.463611.14030-8.00000]}^T,{\mathbf{x}}^{*,2}={[11.070606.79390-2.48650]}^T,\\ {\mathbf{x}}^{*,3}={\left[9.198102.0589015.87207\right]}^T,{\mathbf{x}}^{*,4}={[-0.278002.51633-1.99493]}^T,\\ {\mathbf{x}}^{*,5}={[-3.984302.2704813.40976]}^T,{\mathbf{x}}^{*,6}={\left[17.182350.262532.19246\right]}^T,\\ {\mathbf{x}}^{*,7}={[-3.567843.4484214.94369]}^T,{\mathbf{x}}^{*,8}={[4.538812.20612-2.99328]}^T,\\ {\mathbf{x}}^{*,9}={\left[9.064192.4992811.30274\right]}^T,{\mathbf{x}}^{*,10}={[-0.16333-1.88939-0.75766]}^T.\\ {\mathbf{x}}^{*,11}={[2.179566.92983-1.06412]}^T,{\mathbf{x}}^{*,12}={[10.24873-0.514655.82840]}^T,\\ {\mathbf{x}}^{*,13}={\left[1.121642.845067.93804\right]}^T,{\mathbf{x}}^{*,14}={[0.106783.23489-1.55778]}^T,\\ {\mathbf{x}}^{*,15}={[-2.543740.997322.82005]}^T,{\mathbf{x}}^{*,16}={\left[7.410066.4963412.02799\right]}^T,\\ {\mathbf{x}}^{*,17}={[-0.675104.03845-0.28527]}^T,{\mathbf{x}}^{*,18}={[-0.180374.629806.89661]}^T.\end{array}$$

(51)

A value of the quality criterion (50) for one initial state $\mathbf{x}\left(0\right)={\left[050000\right]}^T$, is $J_4=15.6090$. In the Figure 4 projections of the optimal trajectory on the horizontal plane $\{x_1;x_3\}$ found by the adaptive synthesized control (51) is presented.

Since the initial state in the problem coincided with the terminal state, in order to force the control object to move along a closed path, mandatory conditions for passing through certain areas were added to the quality criterion. For trajectories that meet the criteria for passing through the specified areas, the value of the quality criterion will not change at $p_3=0$. How it is seen in the Fig.3 and the Fig.4 both trajectories pass through all specified areas.

Let’s check the sensitivity of the obtained solutions to random perturbations of the initial state

$$x_i\left(0\right)=x_i^0+{\beta }_0(2\xi \left(t\right)-1),i=1,\dots ,n=6,$$

(52)

where $\xi \left(t\right)$ is a function generator of random noise, returns random number from interval $(0;1)$ at every call, ${\beta }_0$ is a constant level of noise.

In the Figure 5 and 6 the optimal (in blue) and eight perturbed trajectories (in black) for ${\beta }_0=0.1$ for solutions obtained by synthesized (Figure 5) and adaptive synthesized (Figure 6) control are presented.

For comparison, for a model (35) without stabilization systems (39), the problem of optimal control directly was solved, where control was sought in the form of a piece-wise linear function, taking into account restrictions (36).

$$u_i=\left\{\begin{array}{c}{u}_i^{+},\mathrm{if}{\tilde{u}}_i>u_i^{+}\hfill \\ u_i^{-},\mathrm{if}{\tilde{u}}_i<u_i^{-}\hfill \\ {\tilde{u}}_i,\mathrm{otherwise}\hfill \end{array}\right.,i=1,\dots ,m=4,$$

(53)

where

$${\tilde{u}}_i=(q_{i+jm}-q_{i+(j-1)m})\frac{t-(j-1)\Delta t}{\Delta t}+q_{i+(j-1)m},i=1,2,3,4,$$

(54)

$(j-1)\Delta t\le t<j\Delta t$, $j\in \{1,\dots ,K+1\}$, $q_i$ is a component of desired parameters vector, $i=1,\dots ,m(M+1)$,

$$\mathbf{q}={[q_1\dots q_{m(M+1)}]}^T.$$

(55)

In this work, we set the same time interval $\Delta t=0.8$, therefore from (15) $M=18$, and it is necessary to find $m(M+1)=4·19=76$ parameters, $\mathbf{q}={[q_1\dots q_{76}]}^T.$

To solve the optimal control problem the same hybrid algorithm [28] was used. As a result, the following solution was obtained

$$\begin{array}{ccc}\hfill \mathbf{q}& =& [12.57045-4.58471-2.746172.90422-9.26325-0.10990\hfill \\ \hfill & & -2.6322218.368410.00816-17.351770.001654.28718\hfill \\ \hfill & & -11.814921.88489-8.0120615.8794310.848940.06505\hfill \\ \hfill & & 9.6547519.519032.79860-4.06408-0.8899210.50507\hfill \\ \hfill & & -19.1903017.9024012.5243119.000104.76513-11.97648\hfill \\ \hfill & & 0.000108.854642.923340.142388.609197.83194\hfill \\ \hfill & & 5.74904-8.35383-3.4275712.8767118.5871715.43057\hfill \\ \hfill & & 9.0613712.55621-1.546281.473142.407068.67602\hfill \\ \hfill & & 0.00091-11.91236-19.9406317.0830419.92640-1.33145\hfill \\ \hfill & & -7.7725815.54094-19.93278-17.37121-9.312905.03257\hfill \\ \hfill & & -0.90297-5.220210.626534.21368-2.04314-0.53192\hfill \\ \hfill & & 0.0935314.25213-0.115879.05588-0.0327011.23667\hfill \\ \hfill & & {0.03826-16.780470.1822019.81652]}^T.\hfill \end{array}$$

(56)

In the Figure 7 the projection of the optimal trajectory obtained by the direct method is presented.

In the table 1 values of the quality criterion (44) of ten experiments for perturbed solutions obtained by the synthesized (column Synthesized), the adaptive synthesized (column Adaptive) control, and the direct solution (Direct). In two last strings of the table average values of the functionals and standard deviations for all experiments are presented.

As can be seen from Figure 5 and 6 and the Table 1, the solutions obtained by adaptive synthesized control are less sensitive to perturbations of initial states than the solutions obtained by simple synthesized control or especially by the direct approach.

5. Results

A new method for solving the problem of optimal control in the class of implemented functions, an adaptive synthesized control principle is presented. Unlike synthesized control, the new method takes into account the perturbations of the initial state when solving the optimal control problem. Therefore, the value of quality criterion is calculated as the sum of the quality criterion values for the different initial states. As a result of this approach, a solution is chosen in such a way that for the origin initial state it may not give the best quality criterion value, but in the case of disturbances of the initial state, the quality criterion value changes slightly.

6. Discussion

Obtaining a solution based on replacing the optimal solution is less optimal, but also less sensitive to disturbances, at first glance, seems obvious and can be applied to any method of solving the optimal control problem. However, this is not the case. A direct solution to the optimal control problem results in control in the form of a time function and an open loop control system. Perturbation of the initial conditions for such a system gives large variations in quality criterion values, which cannot be reliably estimated from the average value due to the large variance.

The synthesized control method firstly makes the control object stable relative to some equilibrium point in the state space. This means that the perturbed and unperturbed trajectories at each point in time move towards a stable equilibrium point. The adaptive synthesized control method sets the positions of the equilibrium points so that all disturbed trajectories are located in some tube that does not violate phase constraints whenever possible.

In the future, when using the adaptive synthesized control method, it is necessary to assess the required size of the initial state region and reduce the number of initial state points, since this significantly increases the time for finding the optimal solution.