Элементы теории оптимального управления

Размещено на http://www.allbest.ru/

КАМЧАТСКИЙ ГОСУДАРСТВЕННЫЙ ТЕХНИЧЕСКИЙ УНИВЕРСИТЕТ

Кафедра «Системы управления»

Контрольная работа

«Элементы теории оптимального управления»

Петропавловск-Камчатский 2013

Задача № 1

Найти экстремаль функционала:

Исследовать полученную экстремаль на достаточные условия экстремума. Вычислить значение функционала на найденной экстремали и, для сравнения, на прямой, соединяющей точки M1(x1,y1) и M2(x2,y2). Построить график решения.

В этом примере подынтегральная функция F(x,y,y') является функцией общего вида, поэтому составим уравнение Эйлера и решим его. Затем построим график решения. Попутно исследуем на выполнение достаточных условий экстремума и вычислим значение функционала на экстремали и отрезке прямой M1M2.

Применим для решения задачи MATLAB.

Очистим память. Напечатаем заголовок решаемой задачи. Опишем символические переменные. Для решения уравнения Эйлера используем принятые в MATLAB обозначения производных: Dy для y' и D2y для y''. Аргумент обозначим x, а функцию ? y.

>> clear all % Очистили память Workspace

>> disp('Вариант 24 задача № 1') % Выводим заголовок задачи

>> syms x y Dy D2y % Описали символические переменные

Решаем задание:

Вводим подынтегральную функцию и граничные условия. Печатаем их. Чтобы напечатать символ ', нужно его удвоить. Здесь вы должны поставить свои исходные данные: подынтегральную функцию F и граничные условия x1, y1, x2, y2.

>> F=(Dy.^2-9*y.^2)+(4*y*(sin(3*x)))+(5*x.^2); % Подынтегральная функция

>> % Задаем граничные условия

>> x1=-1;

>> y1=2;

>> x2=0;

>> y2=0;

>> disp('исходные данные:')

>> fprintf(['Подынтегральная функция''F(y'',y,x)=%s\n'],char(F))

>> fprintf('Граничное условие слева: y(%d)=d\n',x1,y1)

>> fprintf('Граничное условие справа: y(%d)=d\n',x2,y2)

Исходные данные:

Подынтегральная функция F(y',y,x) = Dy^2-9*y^2+4*y*sin(3*x)+5*x^2

Граничное условие слева: y(-1) = 2

Граничное условие справа: y(0) = 0

Начинаем вывод дифференциального уравнения Эйлера. Найдём частные производные Fy и Fy'. Напечатаем их.

>> dFdy=diff(F,y); % Вычислили Fy

>> dFdy1=diff(F,Dy); % Вычислили Fy'

>> fprintf('Fy=%s\n',char(dFdy))

>> fprintf('Fy''=%s\n',char(dFdy1))

Fy=-18*y+4*sin(3*x)

Fy'=2*Dy

В уравнение Эйлера входит полная производная dFy'/dx. Вычислим её по обычной формуле дифференцирования сложной функции:

Напечатаем её. Напечатаем также величину Fy'y', необходимую для проверки достаточных условий экстремума по признаку Лежандра.

>> d_dFdy1_dx=diff(`dFdy1,x'); % Fxy'x

>> d_dFdy1_dy=diff(`dFdy1,y'); % Fyy'

>> d_dFdy1_dy1=diff(`dFdy1,Dy'); % Fy'y' - условие Лежандра

>> dFy1dx=d_dFdy1_dx+d_dFdy1_dy*Dy+d_dFdy1_dy1*D2y;

>> fprintf('dFy''/dx=%s\n',char(dFy1dx))

>> disp('Условие Лежандра:')

>> fprintf('Fy''y''=%s\n',char(d_dFdy1_dy1))

dFy'/dx=2*D2y

Условие Лежандра:

Fy'y'=2

функционал эйлер экстремаль

Составим левую часть дифференциального уравнения Эйлера и упростим её. Преобразуем символическую переменную Euler в строку.

>> Euler=simple(dFdy-dFy1dx); % Уравнение Эйлера

>> deqEuler=[char(Euler) '=0']; % В строку добавили =0

>> fprintf('Уравнение Эйлера: \n%s\n',deqEuler)

Уравнение Эйлера:

-18*y+4*sin(3*x)-2*D2y=0

Мы составили уравнение Эйлера, теперь решим его. Команда dsolve позволяет находить как общее решение дифференциального уравнения, так и частное его решение, удовлетворяющее заданным начальным или граничным условиям. В следующих главах при решении других заданий нам нужно будет иметь общее решение уравнения Эйлера. Найдём его.

>> Sol=dsolve(deqEuler,'x'); % Решаем уравнение Эйлера

>> if length(Sol)~=1 % Решений нет или более одного решения

error('Нет решений или более одного решения!');

else

disp('Общее решение уравнения Эйлера:')

fprintf('y(x)=%s\n',char(Sol))

end

Общее решение уравнения Эйлера:

y(x)=sin(3*x)*C2+cos(3*x)*C1-1/3*cos(3*x)*x

Сформируем теперь уравнения для граничных условий. Подставим в найденное аналитическое решение Sol граничные точки x1 и x2, и приравняем их соответственно y1 и y2.

>> SolLeft=subs(Sol,x,x1);

>> SolRight=subs(Sol,x,x2);

>> EqLeft=[char(SolLeft) '=' char(sym(y1))];

>> EqRight=[char(SolRight) '=' char(sym(y2))];

>> disp('Уравнения для граничных условий:')

>> fprintf('%s\n',EqLeft,EqRight)

Уравнения для граничных условий:

-sin(3)*C2+cos(3)*C1+1/3*cos(3)=2

C1=0

Решаем полученную систему конечных уравнений ? находим значения произвольных постоянных C1 и C2. Присваиваем найденные решения символическим константам, полученным при решении дифференциального уравнения. Теперь вычисляем аналитическое решение Sol21. Такое вычисление сводится к тому, что в него будут подставлены найденные значения констант C1 и C2. Печатаем найденное уравнение экстремали.

>> Con=solve(EqLeft,EqRight,'C1,C2'); % Решаем систему

% Присваиваем полученные решения символическим константам С1 и С2

>> C1=Con.C1;

>> C2=Con.C2; %

>> Sol21=vpa(eval(Sol),14); % Подставляем С1 и С2

>> disp('Уравнение экстремали:')

>> fprintf('y(x)=%s\n',char(Sol21))

Уравнение экстремали:

y(x)=-16.510752308619*sin(3.*x)-.33333333333333*cos(3.*x)*x

Вычислим значения функционала на найденной экстремали и на прямой, соединяющей точки M1 и M2. Подставим в подынтегральную функцию F аналитические выражения для этих линий и их производных, а затем проинтегрируем. Напечатаем результаты.

>> Fextr=subs(F,{y,Dy},{Sol21,diff(Sol21,x)});

>> Jextr=eval(int(Fextr,x,x1,x2))

>> ylin=(x-x1)*(y2-y1)/(x2-x1)+y1;

>> Flin=subs(F,{y,Dy},{ylin,diff(ylin,x)});

>> Jlin=eval(int(Flin,x,x1,x2))

Jextr = -114.5724

Jlin = -9.0988

В данном примере условие Лежандра говорит о сильном минимуме, что подтверждается полученным результатом: значение функционала на экстремали меньше, чем на другой допустимой функции.

И, наконец, строим график. Задаём массив аргументов для рисования графика функции и вычисляем значения функции. Рисуем график, подписываем заголовок и координатные оси установленным шрифтом.

>> xpl=linspace(x1,x2); % Задаем массив абсцисс

>> y21=subs(Sol21,x,xpl); % Вычислили ординаты

>> figure % Строим фигуру

>> plot(xpl,y21,'-r') % Рисуем график красной линией

>> set(get(gcf,'CurrentAxes'),'FontName',

'Times New Roman Cyr','FontSize',12)

>> title('\bfЗадание 1') % Заголовок

>> xlabel('\itx') % Метка по оси ОХ

>> ylabel('\ity\rm(\itx\rm)') % Метка по оси OY

В итоге получаем:

Вариант 24 задача № 1

Исходные данные:

Подынтегральная функцияF(y',y,x)=Dy^2-9*y^2+4*y*sin(3*x)+5*x^2

Граничное условие слева: y(-1)=d

Граничное условие слева: y(2)=d

Граничное условие справа: y(0)=d

Граничное условие справа: y(0)=d

Fy=-18*y+4*sin(3*x)

Fy'=2*Dy

dFy'/dx=2*D2y

Условие Лежандра:

Fy'y'=2

Уравнение Эйлера:

-18*y+4*sin(3*x)-2*D2y=0

Общее решение уравнения Эйлера:

y(x)=sin(3*x)*C2+cos(3*x)*C1-1/3*cos(3*x)*x

Уравнения для граничных условий:

-sin(3)*C2+cos(3)*C1+1/3*cos(3)=2

C1=0

Уравнение экстремали:

y(x)=-16.510752308619*sin(3.*x)-.33333333333333*cos(3.*x)*x

Jextr = -114.5724

Jlin = -9.0988

Задача № 2

Найти экстремаль функционала:

Исследовать на выполнение достаточных условий экстремума. Построить график решения. В этом примере подынтегральная функция F(x,y,y') не зависит явно от y. Составим программу для решения этой вариационной задачи. Вначале введём исходные данные. У нас будет первый интеграл уравнения Эйлера, поэтому ни сама функция y, ни её вторая производная y'' нам не нужны, и мы их не описываем. Поставим подынтегральную функцию и граничные условия.

>> clear all

>> disp('Вариант 24 задача № 2')

>> syms x Dy

>> F=(Dy^2)+(2*x*Dy*(log(x)))-(log(x));

>> x1=1;

>> y1=-1;

>> x2=3;

>> y2=2;

>> disp('исходные данные:')

>> fprintf(['Подынтегральная функция',

'F(x,y'')=%s\n'],char(F))

>> fprintf('Граничное условие слева: y(%d)=%d\n',x1,y1)

>> fprintf('Граничное условие справа: y(%d)=%d\n',x2,y2)

Решаем задание.

Исходные данные:

Подынтегральная функция F(x,y') = Dy^2+2*x*Dy*log(x)-log(x)

Граничное условие слева: y(1) = -1. Граничное условие справа: y(3) = 2

Строим первый интеграл и решаем полученное дифференциальное уравнение. Названия констант C1 и C2 используются в команде dsolve, поэтому при составлении 1-го интеграла уравнения Эйлера обозначим константу C.

>> dFdy1=simple(diff(F,Dy));

>> deqEuler=[char(dFdy1) '=C']

>> disp('Первый интеграл уравнения Эйлера:')

>> fprintf('%s\n',deqEuler)

>> Sol=dsolve(deqEuler,'x');

>> if length(Sol)~=1

error('Нет решений или более одного решения!');

else

disp('Общее решение уравнения Эйлера:')

fprintf('y(x)=%s\n',char(Sol))

end

Первый интеграл уравнения Эйлера: 2*Dy+2*x*log(x) = C

Общее решение уравнения Эйлера:

y(x) = -1/2*x^2*log(x)+1/4*x^2+1/2*C*x+C1

В переменной Sol получено общее решение, произвольные постоянные обозначены C и C1. Найдём их. Для этого подставим в Sol граничные точки. Приравняем полученные выражения соответственно y1 и y2. Тем самым мы сформируем систему уравнений.

>> SolLeft=subs(Sol,x,sym(x1));

>> SolRight=subs(Sol,x,sym(x2));

>> EqLeft=[char(SolLeft) '=' char(sym(y1))];

>> EqRight=[char(SolRight) '=' char(sym(y2))];

>> disp('Уравнения для граничных условий:')

>> fprintf('%s\n',EqLeft,EqRight)

Уравнения для граничных условий:

1/4+1/2*C+C1=-1 -9/2*log(3)+9/4+3/2*C+C1=2

Решим полученную систему ? найдём произвольные постоянные C и C1. Подставим их в решение Sol. Ограничим решение 14 знаками. Напечатаем уравнение найденной экстремали.

>> Con=solve(EqLeft,EqRight,'C,C1');

>> C=Con.C;

>> C1=Con.C1;

>> Sol22=vpa(eval(Sol),14);

>> disp('Уравнение экстремали:')

>> fprintf('y(x)=%s\n',char(Sol22))

Уравнение экстремали:

y(x)=-50000000000000*x^2*log(x)+25000000000000*x^2+2.9718776495032*x-4.2218776495032

Дальнейшие действия не отличаются от описанных. Рисуем график и вычисляем Fy'y', которая нужна для проверки достаточных условий экстремума по признаку Лежандра.

>> xpl=linspace(x1,x2);

>> y22=subs(Sol22,x,xpl);

>> figure

>> plot(xpl,y22,'-r')

>> set(get(gcf,'CurrentAxes'),

'FontName','Times New Roman Cyr','FontSize',12)

>> title('\bfЗадание 2')

>> xlabel('\itx')

>> ylabel('\ity\rm(\itx\rm)')

>> Leg=diff(dFdy1,Dy);

>> disp('Достаточное условие Лежандра:')

>> fprintf('Fy''y''=%s\n',char(Leg))

>> Fextr=subs(F,{Dy},{diff(Sol22,x)});

>> Jextr=eval(int(Fextr,x,x1,x2))

>> ylin=(x-x1)*(y2-y1)/(x2-x1)+y1;

>> Flin=subs(F,{Dy},{diff(ylin,x)});

>> Jlin=eval(int(Flin,x,x1,x2))

Достаточное условие Лежандра:

Fy'y'=2 Jextr = 10.1715 Jlin = 12.0354

В итоге получаем:

Вариант 24 задача № 2

исходные данные:

Подынтегральная функцияF(x,y')=Dy^2+2*x*Dy*log(x)-log(x)

Граничное условие слева: y(1)=-1

Граничное условие справа: y(3)=2

deqEuler = 2*Dy+2*x*log(x)=C

Первый интеграл уравнения Эйлера:

2*Dy+2*x*log(x)=C

Общее решение уравнения Эйлера:

y(x)=-1/2*x^2*log(x)+1/4*x^2+1/2*C*x+C1

Уравнения для граничных условий:

1/4+1/2*C+C1=-1

-9/2*log(3)+9/4+3/2*C+C1=2

Уравнение экстремали:

y(x)=-.50000000000000*x^2*log(x)+.25000000000000*x^2+2.9718776495032*x-4.2218776495032

Достаточное условие Лежандра:

Fy'y'=2

Jextr = 10.1715

Jlin = 12.0354

Задача № 3

Решить задачу о брахистохроне, соединяющей точки M1(0,0) и M2(4,4).

Нам необходимо найти значение константы C1 и параметра в конечной точке t2 из решения системы уравнений.

Составим программу для решения этого примера.

Вначале введём исходные данные задачи. Подставим свою правую точку.

>> clear all

>> disp('Вариант 24 задача № 3')

>> x2=500;

>> y2=4;

>> fprintf('Правая точка: y(%d)=%d\n',x2,y2)

Решаем задание.

Правая точка: y(500) = 4

Составляем систему уравнений. Левую часть каждого уравнения мы задаём сразу в виде строки. В правой части переводим числа x2 и y2 в их строковые представления с помощью функции num2str. Решаем полученную систему уравнений аналитически. Печатаем решение.

>> eq1=['C1*(t2-sin(t2))=' num2str(x2)];

>> eq2=['C1*(1-cos(t2))=' num2str(y2)];

>> fprintf('Система уравнений:\n%s\n%s\n',eq1,eq2)

>> Sol=solve(eq1,eq2,'C1,t2');

>> C1=eval(Sol.C1);

>> t2=eval(Sol.t2);

>> disp('Найденное решение:')

>> fprintf('Константа C1=%10.5f\n',C1)

>> fprintf('Параметр t2=%10.5f\n',t2)

>> disp('Уравнение брахистохроны:')

>> fprintf('x(t)=%10.5f(t-sin(t))\n',C1)

>> fprintf('y(t)=%10.5f(1-cos(t))\n',C1)

Система уравнений:

C1*(t2-sin(t2))=500

C1*(1-cos(t2))=4

Найденное решение:

Константа C1 = 26.56555

Параметр t2= 18.29366

Уравнения брахистохроны:

x(t)= 26.56555(t-sin(t))

y(t)= 26.56555(1-cos(t))

Рисуем график полученной брахистохроны. Выбираем начало координат в левом верхнем углу с помощью команды axis. Задаём границы по оси Ox, чтобы график занимал всё место на рисунке. Выравниваем масштабы по осям координат, чтобы брахистохрона выглядела неискажённой. Надписываем заголовок и метки осей.

>> t=linspace(0,t2);

>> x=C1*(t-sin(t));

>> y=C1*(1-cos(t));

>> figure

>> plot(x,y)

>> set(get(gcf,'CurrentAxes'),'FontName','Times New Roman

Cyr','FontSize',12)

>> axis ij

>> xlim([0 x2])

>> da=daspect;

>> da(1:2)=min(da(1:2));

>> daspect(da);

>> title('\bfЗадание 3 - задача о брохистохроне')

>> xlabel('\itx\rm(\itt\rm)')

>> ylabel('\itx\rm(\itt\rm)')

В итоге получаем:

Вариант 24 задача № 3

Правая точка: y(500)=4

Система уравнений:

C1*(t2-sin(t2))=500

C1*(1-cos(t2))=4

Найденное решение:

Константа C1= 26.56555

Параметр t2= 18.29366

Уравнение брахистохроны:

x(t)= 26.56555(t-sin(t))

y(t)= 26.56555(1-cos(t))

Размещено на Allbest.ru