Деменков Н.П. - Вычислительные методы решения задач оптимального управления на основе принципа максимума Понтрягина - 2015 (842911), страница 4
Текст из файла (страница 4)
Отрезок [а, Ь] делится точкамиотрезки, не обязательно равные) и каждой точкеветствиеYi.х1< х2 < ... Хпностным отношениемЛу / Лхкоторой определяются значениянаXi ставится в соотИсходя из начальных значений х1 ипроксимации производной у' = dy / dx(у1 , путем апв каждой точке сетки разстроится система уравнений, изYi .В процессе ее решения сеткаXi может перестраиваться (в частности, сгущаться) ; при решениииспользуется матрица-якобиан дf / дх.Простейшая форма обращения к функцииbvp4c:sol = bvp4c(odefun,bcfun,solinit),гдеodefun функция, вычисляющая вектор правых частей;bcfun - функция, вычисляющая вектор граничных условий, двекомпоненты которого представляют собой выражения, обращающиеся в нуль в точках а и Ь соответственно; аргументами функции22bcfun являются Уа и Уь - векторы решения у в точках а и Ь; solinit- выходная структура функции bvpinit.С помощью solinit задаются начальные значения Xi и Yi:solinit = bvpinit(xinit,yinit);при этом заполняются два поля:solinit.x = xinit и solinit.y = y(i,:).Здесь xinit вектор-строка a=xinit(l )<xinit(2)< ...
<xinit(n) = Ь;yinit - гипотетические значения для y(i), которые могут быть заданы в одной из двух форм: в виде вектора, каждая компонентакоторогоyinit(i)копируется в качестве гипотетического решениядля всех точек сетки, т. е.y(i,;) = yinit(i), или в виде функции вформе улюбая точка отрезка [а, Ь], у -= quess(x), где х -вектор,длина которого равна порядку системы дифференциальных уравнений .
Для каждой точки сеткиx(i)вычисляется вектор гипотетического решенияу(i, :)= quess(x(i)).Структура sol, аналогичная solinit, содержит решение краевойзадачи; кроме полей sol.x и sol.y она имеет поле sol.yp, в которомсодержатся значения производной решения (sol.y)' в точках sol.x.Полная форма обращения к функции bvp4c:sol = bvp4c(odefun, bcfun, solinit, options, Р 1, Р2, ...
),гдеoptions -аргумент, позволяющий задавать различные управляющие параметры;Pl , Р2 и т. д. - дополнительные аргументы длявычисления odefun и bcfun. Управляющие параметры задаются путем обращения к функции bvpset, аргументом которой является последовательность пар вида«< ' параметр', значение >», например,options = bvpset('FJacobian', @FJac).Параметр FJacobian задает в аналитическом виде якобиангдеf-функция, вычисляемая вПараметрBCJacobian[:Jodefun.задает в аналитическом виде два якобиа-на=[::] и [::],где Ьс(уа,Уь)-функция, вычисляемая в bcfun.23Неизвестный параметр или вектор неизвестных параметроввводится с помощью функцииbvpinit,форма обращения к которойимеет видsolinit = bvpinit(xinit, yinit, parametrs).Здесьsolinit - структура, в которой кроме полей solinit.x = xinit иsolinit.y = y(i,:) заполняется еще одно поле: solinit.parametrs = parametrs; parametrs - гипотетическое значение неизвестного параметра (вектора).Пример1.2.Решим задачу, рассмотреннуюиспользуя функциюв примере1.1,bvp4c.Решение.
Составим m-файл, содержащий функцию, котораявычисляет правые части:function dxdt = exampl(t,x)dxdt = [х(2); ...3 *t+x(2)/t];а также m-файл, содержащий функцию, которая задает граничныеусловияfunction res = border(ха,хЬ)res = [xa(l)-2 xb(l)-9];Начальные значения координатыцииtвыберем с помощью функlinspace пяти равноотстоящих точек на интервале [1 , 2]:tinit = linspace(l , 2, 5).В качестве начальных приближений вектора решений в этихточках выберем, например,xinitполе структурыxinit = [О 2].
При таком выборе вектораsolinit.x будет заполнено копиями исходноговектора. Такой выбор начального приближения может оказаться неочень удачнь~м. Альтернативный вариант состоит в задании вспомогательной функции, вычисляющей эти компоненты. Начальное приближение влияет на число итераций. Так что выбор начального вектора или вида вспомогательной функции полностью зависит отинтуиции пользователя .Сформируем структуруsolinit:solinit = bvpinit(tinit, xinit).Для решения задачи выполним командуsol = bvp4c(@exampl, @border, solinit)24Если построить решение краевой задачи в виде графика с помощью командыplot(sol.t, sol.x(l,:)),результат вычислений получится не очень точный: слишком малоточек. Эту трудность можно преодолеть с помощью дополнительной функцииdeval.Функциядержащуюся в структуреdeval,sol,используя информацию, состроит интерполяционный сплайнЭрмита для заданного вектора пробных точекtinit,который можетсодержать сколь угодно точек для обеспечения требуемой точности результата.
Обращение к функцииdeval имеет видxinit = deval(sol, tinit),гдеxinit -значение сплайн-функции в пробных точках.Командойtinit = linspace(l, 2)создадим вектор, содержащий по умолчанию100 пробных точек.1.4. Решение оптимизационных задачПакет расширенияOptimization Toolboxпредназначен длярешения оптимизационных задач и систем нелинейных уравнений. Он поддерживает следующие основные методы оптимизации функций ряда переменных: безусловную оптимизацию нелинейных функций; метод наименьших квадратов и нелинейнойинтерполяции; решение нелинейных уравнений; линейное программирование;квадратичное программирование; условнуюминимизацию нелинейных функций; метод минимакса; метод многокритериальной оптимизации.Этот пакет дает возможность решать задачи минимизациифункций, нахождения решений уравнений, задачи аппроксимации (подгонки кривых под экспериментальные данные).
Различные типы таких задач вместефункциями пакетас применяемыми дляOptimization Toolbox приведеныих решенияв табл.1.1.25ТаблицаЗадачи, решаемые средствами пакетаТип задачи1.1Optimization ToolboxФункцияМатематическая записьМАТLАВЗадачи минимизацииСкалярнаяminf(a),(одномерная)а1 < а< а2fminbndаминимизацияБезусловнаяminf(x)fminunc,fminsearchхминимизация(без ограничений)minfт хЛинейноепрограммированиеКвадратичноепрограммированиеAeq х = beq,XL :::; х:::; xu1- хт Нх + f т х при условиях2Ь, Aeq·x = beq, X L :::; х:::; xu,mшхquadprogminfт х при условияхМинимизация принийlinprogАх:::; Ь,Ах:::;наличии ограниче-при условияхххс(х):::; О,ceq(x) =Ах:::; Ь, Aeq·x = beq, XLfminconО,:::; х:::;x u,mш у при условияхх,уДостижение целиF(x)-wy :::; goal,с(х):::; О, ceq(x) = О,Ах:::; Ь, Aeq·x = beq, XL :::; х :::; x u,min max {F;, ( х)}fgoalattainпри условиях{F;(x)}хМинимаксс(х)Ах:::;:::; О, ceq(x) =Ь, Aeq·x = beq, XLfminimaxО,:::; х:::;xu,minfт х при условияхПолубесконечнаяминимизацияхK(x,w):::;с(х):::; О,Ах :::;26О для всехw,ceq(x) = О,Ь, Aeq·x = beq, XL:::; х:::; xu,fseminfОкончание табл.Тип задачи1.1ФункцияМатематическая записьМАТLАВНахождение решений уравненийЛинейные уравненияС(х) ~d , п уравнений, п переменных\ (операторлевогоделения,slashj(a) = ОНелинейноеfzeroуравнение однойпеременнойНелинейныеF(x) = О, п уравнений,уравнения многихfsolveп переменныхпеременныхЗадачи аппроксимации (подгонки кривых)Линейный методminllcx-dll~, т уравнений, пнаименьшиххлевогоделения,квадратов (МНК)переменныхНеотрицательныйminllcx -dll~ при условииbackslash)хлинейный МНК\ (операторlsqnonnegх~ОЛинейный МНК приminllcx - dll~ при условияххналичии ограниченийАх~ Ь,lsqlinAeq·x = beq,XL ~ х ~ хи,min _!_IIF(x)II~ = 2_}i (х)2х 22 l.I:Нелинейный МНКlsqnonlinпри условииXL ~ х ~хи,2Нелинейнаяmin_!_IIF(x,xdata)- ydatallх 22«подгонка» кривойlsqcurvefitпРИ УСЛОВИИ XL ~ Х ~ X u ,В таблице приняты следующие обозначения: агумент; х, у--скалярный арв общем случае векторные аргументы; j{a),скалярные функции;F(x),с(х),ceq(x), K(x,w) -j{x) -векторные функции;27А,С, Н - матрицы; Ь,Aeq,XL, xu -beq, d,f, w, goal, xdata, ydata -векторы;нижняя и верхняя границы области изменения аргумента.Если при поиске минимума вещественной функции у= f ( х)векторного аргумента х на аргумент наложены те или иные ограничения в виде уравнений и/или неравенств, то основным методомрешения задач условной минимизации является использованиемножителей Лагранжа.
Каждое неравенство видаg(x) ~ О превращается в уравнение g(x) + v =О, в котором слагаемое v 2 заведо2мо не отрицательно. Затем левая часть каждого уравнения добавляется к целевой функции с некоторыми множителями, которыеv включаются в число переменных.fgoalattain служит для решения задачи многомернойвместе со значениямиФункциявекторной оптимизации методом «достижения цели» и (при определенных условиях) записывается в следующем виде:•х= fgoalattain(fun,xO,goal,weight) -целевых функцийfun,при заданных вектореначальном приближении хО, заданных векторах целейgoal и весов weight;• х = fgoalattain(fun,xO,goal,weight,A,b) - при наличии ограничений в форме линейных неравенств Ах < Ь;• х = fgoalattain(fun,xO,goal,weight,A,b,Aeq,beq) - при наличии дополнительных ограничений в форме равенств Aeq х = beq;если ограничения в формеА=[]•иЬхнеравенств отсутствуют, задаются= [ ];= fgoalattain(fun,xO,goal,weight,A,b,Aeq,beq,lb,ub) -при наличии дополнительных граничных ограниченийlb ~ х ~ иЬ;• х = fgoalattain(fun,xO,goal,weight,A,b,Aeq,beq,lb,ub, nonlcon) -при наличии дополнительных ограничений в форме нелинейныхнеравенств или равенств с(х) ~ О,ceq(x) = О, если граничные ограничения отсутствуют, то задаются lb = [] и ub = [ ];• х = fgoalattain(fun,xO,goal,weight,A,b,Aeq,beq, ..
.lb,ub, nonlcon,options) - при заданных (изменениях) опциях;• х = fgoalattain(fun,xO,goal,weight,A,b,Aeq,beq, ...lb,ub,nonlcon, options,Pl,P2, ...) - при заданных параметрах Pl, Р2, ... , относящихсяк функциям-аргументам;• [x,fval] = fgoalattain( ...) -возвращает не только оптимальноезначение векторного аргумента, но и значение целевой функции вточке минимума fval;28• [x,fval,attainfactor] = fgoalattain( ...) -то же, что и предыдущая функция, но возвращает еще и коэффициент достижения целиattainfactor;• [x,fval,attainfactor,exitflag] = fgoalattain( ... ) -то же, что ипредыдущая функция, но возвращает еще информацию о характере завершения вычисленийexitflag;• [x,fval,attainfactor,exitflag,output] = fgoalattain( ...