Научная статья на тему 'Symbolic regression based solution for the optimal control problem with constraints'

Symbolic regression based solution for the optimal control problem with constraints Текст научной статьи по специальности «Математика»

CC BY
108
31
i Надоели баннеры? Вы всегда можете отключить рекламу.
Область наук
Ключевые слова
OPTIMAL CONTROL / SYMBOLIC REGRESSION / NETWORK OPERATOR / PHASE CONSTRAINTS / CONTROL SYNTHESIS

Аннотация научной статьи по математике, автор научной работы — Diveev A.I., Hussein O., Shmalko E.Yu.

This article presents a numerical solution to the problem of optimal control of objects in an environment with phase constraints. The proposed approach of synthesized optimal control consists of two steps. First, the problem of synthesizing the stabilization system of an object relative to some point in the state space is solved. The resulting feedback control system is added to the mathematical model of the control object and then the problem of the optimal location of stabilization points, which are essentially attractors, is solved. To solve the synthesis problem, methods of symbolic regression are used, which are completely machine treated, universal and independent of the type of control object. An example of solving the optimal control problem for a group of quadrocopters moving a cargo on flexible rods in a space with constraints is given.

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Текст научной работы на тему «Symbolic regression based solution for the optimal control problem with constraints»

Symbolic Regression based Solution for the Optimal Control Problem with Constraints

Askhat I. Diveev, Oubai Hussein, and Elizaveta Yu. Shmalko

Abstract—This article presents a numerical solution to the problem of optimal control of objects in an environment with phase constraints. The proposed approach of synthesized optimal control consists of two steps. First, the problem of synthesizing the stabilization system of an object relative to some point in the state space is solved. The resulting feedback control system is added to the mathematical model of the control object and then the problem of the optimal location of stabilization points, which are essentially attractors, is solved. To solve the synthesis problem, methods of symbolic regression are used, which are completely machine treated, universal and independent of the type of control object. An example of solving the optimal control problem for a group of quadrocopters moving a cargo on flexible rods in a space with constraints is given.

Keywords—Optimal control, symbolic regression, network operator, phase constraints, control synthesis.

I. INTRODUCTION

The problem of optimal control has a rather long history and today a large pool of both engineering and analytical methods for its solution has been created. But there is a gap between these approaches.

Analytical methods, well formalized and very beloved by fundamentalist mathematicians, are practically not used by practical engineers, since in real problems with nonlinear objects and functionals, with phase constraints, the use of analytical methods is not always a solvable problem. For example, there is a fact that the Pontryagin maximum principle [1, 2], being the main analytical result on optimal control, is practically not used in complex engineering problems. And the reasons for this are different. Firstly, in view of computational complexity, and secondly, this approach does not say anything about the stability of the movement of the object relative to the resulting optimal trajectory. And the solution of the additional problem of stabilization of motion near the optimal trajectory casts doubt on the optimality of the initially obtained optimal trajectory [3]. Existing numerical methods [4, 5] also do not

Manuscript received July 4, 2020. This work was supported in part by the Russian Foundation for Basic Research, Grant 18-29-03061-mk.

A.I. Diveev is with the Federal Research Center "Computer Science and Control" of the Russian Academy of Sciences, Moscow, Russia (e-mail: aidiveev@mail.ru).

O. Hussein is with RUDN University, Moscow, Russia (e-mail: eng .nano.oubai@gmail .com).

E.Yu. Shmalko is with the Federal Research Center "Computer Science and Control" of the Russian Academy of Sciences, Moscow, Russia (corresponding author, phone +7-964-636-6669, e-mail: e.shmalko@gmail .com).

guarantee stable motion with respect to the optimal trajectory and require additional refinement.

Unlike mathematicians, engineers first make the object stable, and only then calculate the optimal trajectory. The introduction of additional regulators in the feedback changes the model of the control object. And this is an essential fact that mathematicians do not take into account when calculating optimal controls.

However, the calculation of optimal regulators by existing methods both analytical [6, 7] or technical [8, 9] is highly dependent on a specific object and is not a universal approach.

In this paper, it is proposed to use a new two-stage computational method of synthesized optimal control. According to this method, optimal control is calculated for a stabilized object, which meets the requirements of engineers. Moreover, the calculation of the stabilization system and further calculation of the optimal trajectory is carried out by numerical methods that are not dependent on the specific model of the object. This approach is universal and allows you to automate the process of calculating control systems and thereby speed up their development.

The article presents a mathematical model of the control object, which is a group of quadrocopters that together on flexible rods carry the load to a given terminal position in an environment with obstacles. The paper provides a mathematical description of the proposed two-stage approach, data for a computational experiment and the results obtained.

II. Synthesized Optimal Control

Consider the classical problem of optimal control

x = f (x, u), (1)

where x e rn, u e U c rm , U is a compact set.

x(0) = x0 e IT. (2)

x(tf) = xf e IT. (3)

where t f < t+, t+ is a positive value, it is given limited

time of control process.

It is necessary to find the control u e U to minimize the functionality

tf

J = f f 0(x(t), u(t))dt ^ min n u(-)eU

(4)

To solve the optimal control problem (1) - (4) by the synthesized optimal control, initially the synthesis problem is solved.

u = h(x - x) e U

(5)

where Vx e X ç Wn 3x eKn

f (x, h(x - x)) = 0, (6)

det(XE - A(x)) = X" + pn-1 (x)Xn + ... + p1 (x)X + p0(x) =

n(Xj(x)-X) = 0 ,

j=1

where

A(x) =

df(x, h(x - x)) dx

(7)

(8)

X j (~) = a j (~) + z'P j (~), (9)

a j (x) < 0, j = 1,..., n, i = V-1. (10)

The equations (5) - (10) show, that after solution of the synthesis problem the differential equation system

x = f (x, h(x - x)),

(11)

tf

J = f f0 (x(t), h(x* (t) - x(t)))dt ^ min.

J *

0 x O

(12)

(14)

problem (2), (3), (11), (12) to solve by a numerical method. Here Pontriagin's maximum principle can be used.

Hamiltonian for the problem (2), (3), (11), (12) is

* *

H (x, x , y) = — fo (x, h(x - x)) +

f y if (x, h(x * — x)), (16)

i=1

where y = [yx ... yn] is a co-state variable vector

*

dH(x, x , y) . ,

yy i =--L , i = 1,., n . (17)

dxi

*

Since x has no restrictions, the conditions of the maximum principle for the Hamiltonian can be obtained analytically from equation

*

dH (x, x , y)

= 0.

is stable according to Lyapunov's stability theorem by the

first approximation. For these system (11) there is a stable

~ *

point x(x ) , that has an attractor property and is a trivial solution of the system (11).

In the second stage the optimal control problem is solved with the model (11) of control object, the initial conditions

(2), the terminal conditions (3), and the initial functional (4),

*

but only it is necessary to find a vector function x (t) that is used instead of a control vector

dx

o *

5 2 H (x, x , y) (dx*)2

> 0,

(18)

(19)

Note, the control vector x (t) in the new optimal control problem (2), (3), (11), (12) has the same dimension as a state space vector x(t) and doesn't have restriction.

The absence of restrictions on control allows us to consider this problem (2), (3), (11), (12) as Lagrange's task

in the classical variational calculus. If the equation (11) is

*

resolved with respect to x

*

x = g(x, x), (13)

and insert this equation in the functionality (12), then the classical functionality of variational calculus is received

The synthesizing function h(x - x) can very complex and therefore Hamiltonian (16) can have some maximums.

For this optimal control problem (2), (3), (11), (12) better to use direct approach and numerical methods of nonlinear programming.

Our experience shows that the best of all to look for the solution of this optimal control problem (2), (3), (11), (12) in the form piecewise constant functions. As a result coordinates of some points

x*,!,...,x*,^ , (20)

are found. These points are switched on a given time interval

x*(t) = x*,j+1 , if jAt < t < (j + 1)At, j = 0,...,K -1, (21) where

K =

r

At

(22)

tf

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

J = j fo (x(t), h(g(x(t), x(t)) — x(t)))dt.

0

For this functionality (14) Euler function can be apply dfo (x(t ), h(g(x(t ), x (t )) — x(t))) — dx

d_ f dfo (x(t), h(g(x(t), x(t)) — x(t))) ^ o

dt [ dx J . )

As a result n differential of second orders are obtained. Solution of these equation for initial conditions (2) and terminal conditions (3) allow to receive optimal trajectories in the state space.

*

In practice resolving nonlinear (11) respect x (t) is a very difficult problem. More real the optimal control

where At is a given time interval.

III. The Synthesis Problem and its Solutions

In the proposed approach the synthesis problem is a most difficult part. An essential feature of this problem is that its solution is a multidimensional function of a vector argument.

u = h(x* - x): r-r* (23)

Now, there are known many approaches to the solution of the synthesis problem. All these methods can be divided into three classes.

A) Analytical methods. These are various methods for solving the Bellman equation. In this case the problem of stability providing (6) - (10) is needed to formulate as an optimal problem. The backstepping integrator [10] and analytical construction of aggregated controllers [11] also belong to analytical methods. These can be other special methods for determined mathematical model of control object, for example, analytical construction of optimal controllers [8] for linear systems with quadratic quality criterion.

B) Technical methods. These methods for synthesis problem solution consider an approach, when the control function (5) is given with accuracy to values of parameters. These methods also include methods for creating control systems based on various controllers, including PI and PID controllers. Note, that using of PID or PI controller increases the dimension of the control object mathematical model. This class of methods also include the neural networks. A neural network is a complex function with known structure but with big number of unknown parameters. Finding of these parameters' values is called learning of neural network.

C) Numerical methods. These methods look for a structure and parameters values of control function (5). Today, all these numerical methods for control synthesis problem can be designed on a base of symbolic regression methods [14]. The founder of these methods is John Koza from Stanford University [15]. All methods of symbolic regression code the mathematical expression in the form of special code, and then they search for the optimal mathematical expression by some special evolutionary genetic algorithm. Application of the genetic algorithm for searching for mathematical expression requires development of special operations of crossover and mutation, therefore all methods use various genetic algorithms.

In this work for the synthesis problem solution the network operator method is used. This symbolic regression method codes a mathematical expression in the form of directed graph of network operator. In computer memory the network operator graph is stored in the form of integer matrix.

For search of solution the network operator method uses the variation genetic algorithm. This algorithm is built on the principle of small variations of basic solution. According to this principle one basic solution is coded. Then small variations of the code are determined. For an integer matrix, small variation is a change of one element. All small variations are codded too. So one possible solution is one set of small variation codes. Genetic algorithm operations are performed on the sets of small variations.

Consider an example of network operator code. Let a mathematical expression be

u = X1 + exp(-^1 X1) sin(q 2 X2 + q\), (24)

where x1, x2, q1, q 2 are arguments of the mathematical expression.

To code this mathematical expression, the following functions are enough:

- functions with one arguments

F = ( fl (z) = z, f2 (z) = - z, f3 ( z) = sin( z),

(25)

f4 (Z) = exP( - functions with two arguments

F2 = (fl(Zl, Z 2) = Zi + Z 2, f2(Zl, z 2) = Zi Z 2) . (26) The graph of the network operator for mathematical expression (24) is presented in the fig. 1. In the graph, arguments of mathematical expression are located in sourece

-nodes. The number of functions with one arguments are located near the arcs. The numbers of functions with two arguments are located in nodes. The numbers of nodes are

located in upper parts of nodes.

9,

Fig. 1. Graph of the network operator for mathematical expression (24) The matrix of the network operator foe the represented graph in the fig. 1 is

"0 0 0 0 1 0 0 0 1" 000020100 000001000 000001000 ¥ = 000020030 , (27)

000002100 000000140 000000021 000000001

To search for a solution in the NOP method, a variational genetic algorithm is used. It applies a principle of small variation of the basic solution [16]. Define small variations of the network operator matrix. A small variation makes small changes in the network operator matrix and results in a new network operator. Define one basic solution in the form of a basic network operator matrix and search for the best solution on a set of small variations of the basic solution.

Any variation can be described by an integer vector of four components

w = [w1 w2 w3 w4 ]T , (28)

where Wj is a number of the small variation, w1 e {0,1,2,3} , w2 is a row number of the NOP matrix, w3 is a column number of the NOP matrix, w2,w3 e{1,L}, w4 is a number of the function with one or two arguments, if w1 = 1,

then

E{1,...,|F2|} else w4 e{1,...,|FJ}.

Each new solution in the neighborhood of the basic solution is described by an ordered set of vectors of variations (28)

W = (w1,..., wl ),

(29)

where l is a size of the neighborhood of the basic solution.

We obtain a new network operator matrix in the neighborhood of the basic solution after small variations of the basic matrix according to vectors from the set (29)

¥ = Wo = w'

V o ¥0

(30)

where is a basic network operator matrix.

Each ordered set (29) of vectors of variations defines one of possible solutions in the neighborhood of the basic solution. The search area depends on the size l of the neighborhood of the basic solution and the possibility to replace from time to time the basic solution on the best solution found by the time.

To make crossover, select two possible solutions

W'' = (wu,..., w'-' ), Wj = (w7-1,..., w1,1 ), (31)

Randomly choose a crossover point r e{1,...,'} and

w

exchange the subsets of the vectors of variations from the crossover point till the end. We obtain two new possible solutions

W1 = (wy,..., w'-r-1, wj ■r,..., wj ),

V-t^i,1 wj>-! wi,r wi,l) (32)

, wy ).

W1 = (w1

To make mutation, select in the possible solution W = (W',..., Wl) a mutation point ¿e (1,...,l} and produce a new vector of variations Wu in the position u .

IV. Computational Example

Consider the problem of optimal control for a group of three quadcopters performing the joint task of transporting some load on flexible rods from a certain initial position to a given terminal one in space with phase constraints. The control object is described by a system of nonlinear differential equations.

^j = ^j +(x5 sin x( + xj cos xj ) sin xlj / xj = (x5 sin x{ + xj cos x{ ) / cos xj, X3j = (xj sin xj + xj cos x{ ),

cos x

x4 = xjxj (12 -/3)/11 + M1 //„

xi = xx ( /3 - I, ) /12 + M22 /12,

xj = xjxi ( 11 - /2 ) / /3 + M3 / /3,

(33)

moments of inertia of the quadcopter around the axes xi , i = 1,2,3.

The position of the load is determined by the state vector y e r of the center of mass of the load in the subspace {x

73 x8, x9 }

y = [y y , (34)

where y1, y2, y3 are the positions of the center of mass of the cargo along the axes x7, x8, x9 respectively. Given the initial and final position of the load

r 0 0 0-, T

y = [ y y 2 y3 ] ,

y = [yf y{ y(]T •

(35)

(36)

We suppose that all three quadrocopters are similar, and the rods by which they move the load have the same length:

d1 = d2 = d3 = d.

For equal load distribution on quadrocopters, all three quadrocopters should always be at the same height and be located at the vertices of an equilateral triangle. The size of this triangle determines the height of the load. Let the size of the triangle be determined by the radius of the circumscribed circle R0. Then the height difference between the

quadrocopter xj, j = 1,2,3 and the load y2 is defined as

xj -y =Vd2 -R02 , j = 1,2,3 . (37)

Projections of the reaction force of the load on each quadrocopter are determined from the equations:

w = - w

9 12 >

xj0 = Fj sin xj cos xj cos x( + sin x( sin x2j - w(,

xj = F1 cos x3j cos xj cos xj - g - wj,

xj = Fj sin x( cos xj - cos xj sin x2j sin xj - w]j,

where j - is the quadrocopter number; wf are components

of the weight load of the quadrocopter j from the weight of

the carried load, k = 1,2,3 ; w/ , w2 , w^ are projections of the load on the axis ^7, xx*8, x respectively; x1j , x3j are angles of rotation of the quadrocopter j around horizontal axes; x2j is an angle of rotation of the quadrocopter j around a vertical axis; xJj and xj are angular velocities of rotation of the quadrocopter j around the horizontal axes; x5j is an angular velocity of rotation of the quadrocopter j around a vertical axis; x7j , x9j are coordinates of the center of mass of the quadrocopter j on the horizontal plane; x8j is the height of the quadcopter j; xj0, xjl, xj2 are projections of the linear speed of the quadcopter j on the corresponding axis; M/ are control moments created by propellers of the quadcopter j around the axes xi , i = 1,2,3 ; FJ is the total thrust of all four quadcopter propellers , taking into account the mass correction * ; g is a constant value of the acceleration of gravity; It are

1R cos(r) / d,

w1 = w' k d2 - r2 / d,

w3 = |w'| R sin(r)/ d,

(38)

2 2 w = w

R0 (cos(r) + V3sin(r))/2d,

w22 = w 2 L/ d2 - R2 / d,

(39)

2 2 w = w

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

w,3 = w3

R0 (sin(r) - -v/3 cos(r)) / 2d, R0(cos(r) - >/3sin(r))/2d,

= w3Jd2 - R2 / d,

w3 = w3

R0 (sin(r) + ^ cos(r)) / 2d,

where w

, V(wi)2 - (w2)2 + (w3)2, i = 1,2,3^

R0 = ^3((x1 -x2)2 + (x9 -x2)2)/3,

x„ (x9 + x9 + x9 ) / 3

t = arctan

(40)

(41)

(42)

(43)

v x1 — (x7 + x72 + x3 ) / 3 The problem has three types of phase constraints. 1) Limitations caused by interaction in the group: quadrocopters must be at the vertices of an equilateral triangle.

X = |(x7 — x72)2 + (x1 — x92)2 — (x — x7)2 — (x9 — x93)2\ — sp < o, X2 = |(x7 — x72)2 + (x9 — x9)2 — (x72 — x73)2 — (x9 — 4)2\ — sp < o,

^3 = I (x! - x§ )2 - (x! - x^ )21 - ^ < 0, (44)

^4 = |( x - x8)2 - (x2 - x8)2| < 0,

where sp is a given small positive number.

2) The distance between each quadcopter and the load cannot be greater and should not be significantly less than the length of the rods d .

(x*) = L -d < 0,

following control function

Mt =

M+, if Mt > M+ M", if M, < M" , i = 1,2,3 M,, otherwise

F =

®3 + i(xi) = d-L-s, <0,

(45)

(46) i = 1,2,3,

F +, if F > F+ F -, if F < F -F, otherwise

(50)

where

where L = ^(x7 - >1)2 + (x8 - y2)2 + (x9 - y^)2 er is a given small positive number,

y = x} L^ / 3 cos(Q),

y2=x:§- ~

y = x - (x1 -^,2 / 3 sin(Q)),

Q = arctan(| x72 - x73 | /1 x92 - x93 |),

L1,2 = (x1 - ) + (x1 - ) . 3) There are obstacles in the area of moving of quadrocopters and load.

, __(47)

= rj-yl(x7- ^,j)2+ (x9- x2,j)2 < 0,

^ IN (x ) = 6 + (, - 1)5 + Jy J

The quality of control is determined by a nonlinear functional that includes speed components and penalties for violating phase constraints.

tf 4+4s

J = tf +aj £ ( y (tf ) - yf )2 + a f £ ^j (x1, x2, x3 ))dt +

0 J=1

«3 f £ W (xJ)) + ^(®3+i (xJ)) + 6 + (i-1)5 + j (xi ))

V'-1 V

J=1

dt ^

//

^ min

uJ (■),J=1,2,3 (48)

where ai are weighting factors, i = 1,2,3, uJ (■) = [MJ (■) M2 (■) M3J (■) FJ (Of , J = 1,2,3; 3(A) is a Heaviside function

S( A) =

1, if A > 0 0 - otherwise

It, if t < t+ h Lc(t) <sf, 11 + - otherwise

where

Lc (t) =JZ ( y, (t )-yf )2.

J j=1

All quadcopters must be located at the same height in the corners of the equilateral triangle. The equilateral triangle can change its size and rotation angle relative to a fixed coordinate system Two parameters a circumradius R0 and an angle of equilateral triangle turn t are a control for group of robots.

For solution of the synthesis problem the network operator method [13] was used. As a result, it has got the

M1 = sgn(x6 - X6 ) log(| q6 (x6 - X6) | +1) +

f f 3 * 3/ *

q4(x4 -X4) + (x4 -X4) + q1(x1 -X1) + 3q1(x1 -x^

3 f * f 3 3 *

yq4(x4 -x4)+q1(x1 -x1)+(x4 -x4) +vq1(x1 -x1) +

q2(x2 - x2) + q1(x1* - x1) + s8n(x4 - x4^q4(x4 - x4)

+ (x{ - X4)3,

2 = sgn(sgn(A1)exp(| A1 |) -1) X

Vl sgn(A1)exp(|A1|) -1|, M3 = sgn(x5f - X5 ) log(| q5 (x{ - X5 ) | +1) + -q3(x3 -x3) + x3 + x3 + + q2(x* -x2) +

tanh(0.5B^ + q6 (x6f - x6),

A1 = q3(x3 - x3) + x3 - x3 +

sgn(x5 - x5)log(| q5(x5 - x5) | +1) ,

B1 = q 2( x *- x 2)+q6(xf - x6) + (q3(x*- x3) +

* f f 3

x3 - x3 + sgn(x5 - x5 ) log(| q5 (x5 - x5) | +1)) ,

x, =

i , if x, > x, i, if < x~

i = 1,2,3

~i , otherwise

~1 = A2 q9(x9f - x9 ) cos(xn ) exp(-q12 ) x ^A2log(| q9(xf - x9)cos(xn) |),

~2 = + 2 arctan(C 2) - qf0 x190 + D2 2ff

q11x11q8(x8 - x8) - q7(x7 - x7) +

arctan(A2 q9 (x^ - x9 ) cos(x11 ) exp(-q12 )) + arctan(-q10 xj^ + q7 (x{ - x7 )) + q8 (x8f - x8) +

sgn(-q10x30 + (x/ -+ tanh(0.5x12 ) + sgn(B2 )(exp(| B2 |) -1) + KA2 q9 (x9f - x9 ) cos(xn ) exp(-q12 )),

-H H1 H1 I-

~3 = sgn(~1 ) log(| ~1 | +1) + 3 B2 + arctan(C 2 ) -

3 9 3 f 3

q10x10 - x1 + (-q10x10 + q7(x7 - x7)) , F = sin(~3 ) + sgn(~2* )(exp(| | -1)) + C2 + sgn(~1 ) log(| 5?! | +1) + 2 tanh(0.5B2 ) +

1

f 2 (A2 q9(x9 - x9)cos(x11)exp(-q12)) +

f 3 2

(q7(x7 - x7) - q10 xl0) ,

A2 = q12x12 + + arctan(q9) + cos(x7f -xn),

B2 = arctan(-q10 x^ + qy (x7f - x7 )) -

f 2 f 2 f

q7(x7 - x7) - q11 x11q8(x8 - x8) + q8(x8 - x8), C2 = sgn(-q10 x^ + q7(x7f - x7 )) x (exp(| -q10 xfo + q7(x7f - x7) |) -1) + q7 +

sgn(-q10 x130 ^ -q10 x1301 + ^(q7(x7 - x7)),

D2 = 2(arctan(-q10 x^ + q7(x^ - x7)) -

2 f 2 f q11 x11q8(x8 - x8) + q8(x8 - x8) -1 +

sgn(-qu xn q8 (xf - x8)2) x

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

2 f 2 exp(-1 -qux11q8(x8 -x8) |)+exp(q10)-

q7(x{ -xj) + (-q10x130 + q7(x7f -x7))3) ,

tanh(a) = 1 -exp(2a) , ^(P) = sgn(P) min{| p |,1}, 1 + exp(2a)

q1 = 12.224, q 2 = 14.197, q3 = 13.611, q 4 = 4.361, q5 = 9.989, q 6 = 4.144, q7 = 0.115, q8 = 3.371, q9 = 3.076, q10 = 0.144, q11 = 3.131, q12 = 4.515. On the second stage the vectors of parameters were found

(51)

Q = {q1,., q K }

where

1 f f f T

q = [x7 x8 x9 q4 q5 q6 q7] ,

(52)

x7f , xf , x^ determinate position of the first quadcopter in the geometric space, qi , i = 4,5,6,7 , define location of the equilateral triangle or location of others quadcopters

R =

R0 - R0 ^

. , 2nt 1 . f 4nt , „ sinl--+ q 4 1 + sinl--+ q5 1 + 2

(53)

+ - f T -T

T = -

4

. , 2nt 1 . f 4nt . „ sinj—+ q 6 1 + sinl—+ q7 1 + 2

. (54)

the on-board computer of a real object. Based on the foregoing, we think this approach is very promising.

' x9

Fig. 2 shows the result of simulation. In the fig. 2 black lines are moving trajectories of quadcopters, dot-lines are trajectory of the load, red circles are obstacles.

V. Conclusion

In the paper we described and demonstrated the application of the new two-step approach to the numerical solution of the optimal control problem with phase constraints based on symbolic regression. This approach has several advantages over analytical and technical methods. Firstly, it does not depend on the characteristics of the model; it requires only a mathematical description of the model. The obtained control does not need to be further developed, as the analytical approach requires the stabilization system to be built up. The resulting solution is machine-made and easily transferred to

-7 -6 -5-4 -3 -2 -1 0 1 2 3 4 5 6 7 8 9 Fig. 2. Projections of movement trajectories of three quadcopters and a load on the horizontal plane

REFERENCES

[1] R. Chertovskih, N.T. Khalil, F.L. Pereira, D. Karamzin "Path-constrained trajectory time-optimization in a three-dimensional steady flow field" // 18th European Control Conference, ECC 2019, pp. 3746-3751.

[2] A.V. Arutyunov, D.Yu. Karamzin "Properties of extremals in optimal control problems with state constraints" // Differential Equations, 52(11), 2016, pp. 1411-1422

[3] L.C. Young Lecture on the calculus of variations and optimal control theory, AMS Chelsea publishing, Providence Rhode Island, 2000.

[4] N.K. Khalil "Nonlinear Systems", 3rd Ed., published by Pearson Education, Inc. publishing as Prentice Hall, 1996, 750 p., 2002.

[5] Rao A.V. "A Survey of Numerical Methods for Optimal Control" //

Advances in the Astronautical Sciences, January 2010, 135 (1), pp. 1-32.

[6] Krstic M., Kanellakopoulos M., Kokotovic P.V. Adaptive nonlinear control without overparametrization // Systems & Control Letters. 1992. Vol. 19, Iss. 3. P. 177—185.

[7] A.D. Mizhidon On a Problem of Analytic Design of an Optimal Controller // Automation and Remote Control, 2011, vol. 72, no. 11, pp. 2315-2327.

[8] M. Cutler, J. P. How "Actuator Constrained Trajectory Generation and Control for Variable-Pitch Quadrotors" // A1AA Guidance, Navigation, and Control Conference (GNC). Minneapolis, Minnesota, August 2012.

[9] S.L. Zenkevich, N.K. Galustyan "Synthesis and testing of an algorithm for controlling the motion of a quadrocopter along a trajectory. // Mekhatronika, avtomatizatsiya, upravleniye. 2015; 16(8): 530-535. (in Russian)

[10] M. Krstic, I. Kanellakopoulos, P.V. Kokotovic "Nonlinear and Adaptive Control Design", Wiley, New York, 1995

[11] Podvalny S.L., Vasiljev E.M. Analytical synthesis of aggregated regulators for unmanned aerial vehicles // Journal of Mathematical Sciences, vol. 239, no. 2, May, 2019. pp. 135-145.

[12] Alibekov E., Kubalik J., Babuska R. "Symbolic Method for Deriving Policy in Reinforcement Learning", Proc. 2016 IEEE 55th Conference on Decision and Control (CDC), December 12-14, Las Vegas, USA, 2016, pp 2789-2799.

[13] A.I. Diveev, E.Yu. Shmlako, Self-adjusting control for multi robot team by the network operator method, 2015 European Control Conference (ECC), July 15-17, 2015. Linz, Austria

[14] A.I. Diveev "Numerical methods for solving the control synthesis problem", RUDN (in Russian), 2019.

[15] J.R. Koza. Genetic Programming: On the Programming of Computers by Means of Natural Selection. Cambridge, Massachusetts, London, MA: MIT Press, 1992. -- 819 p.

[16] A.I. Diveev "Small Variations of Basic Solution Method for Non-numerical Optimization", IFAC-PapersOnLine, 2015, volume 48(25) pp. 028-033.

4

t

t

t

t

Askhat I. Diveev was born in Magadan, USSR in 1954. He received the Ph.D. degree in technical science from Bauman Moscow State Technical University, in 1989, and Doctor of Sciences in 2001 in Dorodnitsyn Computing Center of the Russian Academy of Sciences, in 2009 he became a Professor.

He works as a Head of the Department of Robotics Control in the Computing Center of the Federal Research Center "Computer Science and Control" of the Russian Academy of Sciences. He is also a Professor at the RUDN University, Engineering Department. He is the author of five books, more than 300 articles, and seven numerical methods in evolutionary computation. His current research interests are computational methods in control, symbolic regression and evolutionary computation with applications to model identification, optimization, control system synthesis and machine learning.

Prof. Diveev is a member of the editorial board of the RUDN journal of Engineering Researches and journal of Instrument Engineering of the Bauman Moscow State Technical University, a general chair of the INTELS Symposium.

Hussein Oubai was born in damask, Syria, in 1991. She received the B.S. and M.S. degrees in mechanical engineering from the Damask Technical University, Engineering Dept.

He is a graduate student, Department of Mechanics and Mechatronics, Engineering Academy, RUDN University.

Elizaveta Yu. Shmalko was born in Moscow, Russia, in 1982. She received the B.S. and M.S. degrees in computer science and cybernetics from the Russian Peoples Friendship University, Engineering Dept. and the Ph.D. degree from Dorodnicyn Computing Center of Russian Academy of Sciences, Moscow, Russia, in 2009.

She was with IBM East Europe/Asia from 2007 to 2010. Since 2010, she is a Senior researcher with the Computing Center of the Federal Research Center "Computer Science and Control" of the Russian Academy of Sciences, the Department of Robotics Control. From 2013 to 2016 she is also an Assistant Professor at the RUDN University, Engineering Department. She is the author of more than 50 articles, and five evolutionary methods. Her research interests include evolutionary algorithms and symbolic regression methods in application to modeling and control system design in robotic systems and machine learning.

i Надоели баннеры? Вы всегда можете отключить рекламу.