Диссертация (1144294), страница 15
Текст из файла (страница 15)
. . . . . . . . . . . . . . . . . . . . . . . . . . . .1.2Последовательность режимов системы (1.21) и соответствующиеточки переключения. . . . . . . . . . . . . . . . . . . . . . . . . . .1.331Координаты точки на периодическом решении системы (1.21) при = 0.03 и продолжительность соответствующих режимов I и II. .2.13033Значения параметров численной процедуры сканирования области (, ) ∈ (0, 1.1] × (0, 2 + ). . . . .
. . . . . . . . . . . . . . . .67100Приложение АОпределения:динамическиесистемыснепрерывным временем и аттракторыСледуя работам [104,154,159], введем понятие динамической системы заданной автономной системой дифференциальных уравненийx= (x),(А.1)где : R → R — непрерывная вектор-функция, удовлетворяющая локально-му условию Липшица в R . Т. о., по теореме Пикара (см., например, [160, 161])для любого x0 ∈ R существует единственное решение x(,x0 ) дифференциаль-ного уравнения (А.1) с начальным условием x(0,x0 ) = x0 , которое задано нанекотором конечном интервале: ∈ ⊂ R. По теореме о непрерывной зависимости от начальных данных x(,x0 ) непрерывно зависит от x0 [160, 161].При изучении предельного поведения системы (А.1) решения рассматриваются при → +∞ или при → ±∞. В общем случае, для произвольныхквадратичных систем существование решений при ∈ [0 , + ∞), не влечетсуществование при ∈ (−∞, 0 ].
Известно [162], что если непрерывно диф-ференцируемая ( ∈ 1 ), то является локально Липшицевой в R . Такжеизвестно [163], что если функция : R → R является локально Липшицевой,тогда для любого x0 ∈ R решение x(·, x0 ) : → R уравнения (А.1) существуетна максимальном временном интервале = (− , + ) ∈ R, где −∞ ≤ − < 0 и0 < + ≤ +∞. Если + < +∞, тогда |x(, x0 )| → ∞ при → + и, если − > −∞,101тогда |x(, x0 )| → ∞ при → − . Т.о., что решение дифференциального уравне-ния (А.1) продолжимо до тех пор, пока оно остается ограниченным.Для удобства введем множество значений времени T ∈ {R,R+ }. Существо-вание и единственность решений (А.1) при всех ∈ T обеспечивается, например,глобальной Липшецевостью вектор-функции в R .
Еще одним эффективнымметодом для доказательства ограниченности решений при ∈ T является по-строение функций Ляпунова.Известно [160,161,164], что, если выполнены условия существования и единственности для всех ∈ T, то для решения (А.1) справедливо следующее груп-повое свойствоx( + , x0 ) = x(, x(,x0 )),∀ , ∈ T(А.2)и по теореме о непрерывной зависимости решения дифференциального уравнения от начальных данных x(·, ·) : T × R → R является непрерывным отображением.
Таким образом, если решения системы дифференциальных уравне-ний (А.1) существуют и удовлетворяют групповому свойству (А.2) при всех ∈ T, то система (А.1) порождает динамическую систему [165] на фазовом√︀21 + · · · + 2 евклидова норма векторапространстве (R , | · |). Здесь |x| =x = (1 , . . . , ) ∈ R , порождающая метрику в R . Далее для сокращениявместо динамическая система, порожденная дифференциальным уравнением(А.1) будем писать динамическая система (А.1). В качестве начального момента времени естественно выбиратьx(0, x0 ) = x0 .В теории динамических систем при изучении предельного поведения траекторий возникает понятие аттрактора. Далее дадим определение аттрактора,следуя работам [104, 105, 166–169].Определение 9.
Будем говорить, что множество ⊂ R, ̸= ∅ являетсяинвариантным для динамической системы (А.1), если x(, ) = ,Здесь x(, ) = {x(, x0 ) | x0 ∈ } .∀ ≥ 0.Определение 10. Будем говорить, что для динамической системы (А.1) инвариантное, замкнутое и ограниченное множество является:1021. локальным аттрактором, если это минимальное локально притягивающее множество (т.е. lim→+∞ dist(, x(,x0 )) = 0 ∀x0 ∈ (), где ()– некоторая -окрестность множества );2.
глобальным аттрактором, если это минимальное глобально притягивающее множество (lim→+∞ dist(, x(,x0 )) = 0 ∀x0 ∈ R ).Здесь dist(,x) = inf y∈ ||y − x|| обозначает расстояние от точки x ∈ R домножества .Замечание 1. Иногда в определении аттрактора свойство замкнутостиопускается (см. [170]).
Кроме того, также иногда опускается свойство ограниченности (см. [171]). Неограниченные аттракторы возникают при изучении неавтономных систем дифференциальных уравнений в расширенном фазовом пространстве. Отметим, что если динамическая система определенадля ∈ R, то локально притягивающее инвариантное множество состо-ит только из целых траекторий, т.е. если x0 ∈ , то x(, x0 ) ∈ для ∀ ∈ R(см. [169]).Определение 11.
Для аттрактора бассейном притяжения называетсямножество начальных состояний ℬ() ⊂ R , такое чтоlim dist(, x(, x0 )) = 0,→+∞∀ x0 ∈ ℬ().Важную роль в доказательстве существования глобального аттрактора вдинамической системе играет свойство диссипативности. Оценки области диссипативности дают аналитическую локализацию аттрактора в фазовом пространстве.
Таким образом, диссипативность системы с одной стороны, доказывает тот факт, что в фазовом пространстве системы нет траекторий, уходящихна бесконечность при → +∞, и, с другой стороны, дает возможность ука-зать конкретные границы области, в которую входят все траектории системыи оценить скорость вхождения.Определение 12. Множество 0 ⊂ R называется поглощающим для динамической системы (А.1), если для любого x0 ∈ R существует такое0 = 0 (x0 ), что x(, x0 ) ∈ 0 для любого ≥ 0 .103Определение 13 ( [169]). Динамическая система(А.1) называется диссипа-тивной, если она обладает ограниченным поглощающим множеством.Замечание 2. В качестве поглощающего множества в R можно рассматривать [172] шары = {x ∈ R : |x| < }. В этом случае, если существуеттакое > 0, чтоlim sup |x(, x0 )| < ∀ x0 ∈ R ,→∞то говорят, что динамическая система диссипативна в смысле Левинсона.Число называется радиусом диссипативности.Теорема 8 ( [169]).
Если динамическая система(А.1) является диссипатив-ной, то она обладает глобальным аттрактором.Эффективным методом доказательства диссипативности является построение специальных функций Ляпунова [173, 174].Изучение аттракторов автономных систем обычно начинается с анализа состояний равновесия, которые всегда можно найти численно или аналитически.При этом, для локализации аттрактора важно определить его бассейн притяжения.
Таким образом, с вычислительной точки зрения, естественно предложитьследующую классификацию аттракторов, основанную на сложности определения в фазовом пространстве бассейна притяженияОпределение 14 ( [81,175–177]). Аттрактор будем называть самовозбуждающимся, если бассейн притяжения этого аттрактора пересекается с любы-ми открытыми окрестностями его стационарных точек. В противном случае будем называть аттрактор скрытым.104Приложение БРеализация алгоритма построения контрпримера к гипотезе КалманаБ.1ЛокализацияпериодическогорешениявсистемеФиттса с разрывной правой частью методом точечных отображенийЛистинг Б.1:symSolOde1.m – функция, задающая точное решение системы(1.21) в области Σ+ .1function sol = symSolOde1 ( m1 , m2 , b )23syms X_1 X_2 X_3 X_4 real456789101112131415C1 = - ( ( b ^2 + m2 ^2) * X_1 + 2 * b * X_2 + X_3 ...- 1 / ( m1 ^2 + b ^2) ) / ( m1 ^2 - m2 ^2) ;C2 = ( ( m1 ^2 + b ^2) * X_1 + 2 * b * X_2 + X_3 ...- 1 / ( b ^2 + m2 ^2) ) / ( m1 ^2 - m2 ^2) ;C3 = - ( b * ( b ^2 + m2 ^2) * X_1 + ...(3* b ^2 + m2 ^2) * X_2 + 3 * b * X_3 + X_4 - ...b / ( m1 ^2 + b ^2) ) / ( m1 * ( m1 ^2 - m2 ^2) ) ;C4 = ( b * ( m1 ^2 + b ^2) * X_1 + ...( m1 ^2 + 3* b ^2) * X_2 + 3 * b * X_3 + X_4 - ...b / (b ^2 + m2 ^2) ) / ( m2 * ( m1 ^2 - m2 ^2) ) ;syms t positive1617x1 ( t) = 1 / ( ( m1 ^2 + b ^2) * ( b ^2 + m2 ^2) ) + ...105181920212223242526272829303132C1 * exp ( - b * t ) * cos ( m1 * t ) + ...C2 * exp ( - b * t ) * cos ( m2 * t ) + ...C3 * exp ( - b * t) * sin ( m1 * t ) + ...C4 * exp ( - b * t ) * sin ( m2 * t ) ;x2 ( t) = - exp ( - b * t ) * ( ( C1 * b - C3 * m1 ) * cos ( m1 * t ) + ...( C2 * b - C4 * m2 ) * cos ( m2 * t ) + ...( C1 * m1 + C3 * b ) * sin ( m1 * t ) + ...( C2 * m2 + C4 * b ) * sin ( m2 * t ) ) ;x3 ( t) = exp ( - b* t ) * ( ...( -( m1 ^2 - b ^2) * C1 - 2 * m1 * b * C3 ) * cos ( m1 * t ) + ...( ( b ^2 - m2 ^2) * C2 - 2 * b * m2 * C4 ) * cos ( m2 * t ) + ...( 2 * m1 * b * C1 - ( m1 ^2 - b ^2) * C3 ) * sin ( m1 * t ) + ...( 2 * b * m2 * C2 + ( b ^2 - m2 ^2) * C4 ) * sin ( m2 * t ) ) ;x4 ( t) = exp ( - b* t ) * ( ...( (3 * m1 ^2 - b ^2) * b * C1 - ( m1 ^2 - 3 * b ^2) * m1 * C3 ) * cos ( m1 *t ) + ( -( b ^2 - 3 * m2 ^2) * b * C2 + (3 * b ^2 - m2 ^2) * m2 * C4 ) *cos ( m2 * t) + ( ( m1 ^2 - 3 * b ^2) * m1 * C1 + (3 * m1 ^2 - b ^2) * b* C3 ) * sin ( m1 * t ) + ( -(3 * b ^2 - m2 ^2) * m2 * C2 - ( b ^2 - 3 *m2 ^2) * b * C4 ) * sin ( m2 * t ) ) ;3334sol = [ x1 ( t ) ; x2 ( t ) ; x3 ( t ) ; x4 ( t) ];3536endЛистинг Б.2:symSolOde2.m – функция, задающая точное решение системы(1.21) в области Σ+ .1function sol = symSolOde2 ( m1 , m2 , b )23syms X_1 X_2 X_3 X_4 real4567891011121314C1 = - ( ( b ^2 + m2 ^2) * X_1 + 2 * b * X_2 + X_3 + ...1 / ( m1 ^2 + b ^2) ) / ( m1 ^2 - m2 ^2) ;C2 = ( ( m1 ^2 + b ^2) * X_1 + 2 * b * X_2 + X_3 + ...1 / (b ^2 + m2 ^2) ) / ( m1 ^2 - m2 ^2) ;C3 = - ( b * ( b ^2 + m2 ^2) * X_1 + ...(3* b ^2 + m2 ^2) * X_2 + 3 * b * X_3 + X_4 + ...b / ( m1 ^2 + b ^2) ) / ( m1 * ( m1 ^2 - m2 ^2) ) ;C4 = ( b * ( m1 ^2 + b ^2) * X_1 + ...( m1 ^2 + 3* b ^2) * X_2 + 3 * b * X_3 + X_4 + ...b / (b ^2 + m2 ^2) ) / ( m2 * ( m1 ^2 - m2 ^2) ) ;10615syms t positive1617181920212223242526272829303132x1 ( t) = - 1 / ( ( m1 ^2 + b ^2) * ( b ^2 + m2 ^2) ) + ...C1 * exp ( - b * t ) * cos ( m1 * t ) + ...C2 * exp ( - b * t ) * cos ( m2 * t ) + ...C3 * exp ( - b * t) * sin ( m1 * t ) + ...C4 * exp ( - b * t ) * sin ( m2 * t ) ;x2 ( t) = - exp ( - b * t ) * ( ( C1 * b - C3 * m1 ) * cos ( m1 * t ) + ...( C2 * b - C4 * m2 ) * cos ( m2 * t ) + ...( C1 * m1 + C3 * b ) * sin ( m1 * t ) + ...( C2 * m2 + C4 * b ) * sin ( m2 * t ) ) ;x3 ( t) = exp ( - b* t ) * ( ...( -( m1 ^2 - b ^2) * C1 - 2 * m1 * b * C3 ) * cos ( m1 *t ) + ...( ( b ^2 - m2 ^2) * C2 - 2 * b * m2 * C4 ) * cos ( m2 * t ) + ...( 2 * m1 * b * C1 - ( m1 ^2 - b ^2) * C3 ) * sin ( m1 * t ) + ...( 2 * b * m2 * C2 + ( b ^2 - m2 ^2) * C4 ) * sin ( m2 * t ) ) ;x4 ( t) = exp ( - b* t ) * ( ...( (3 * m1 ^2 - b ^2) * b * C1 - ( m1 ^2 - 3 * b ^2) * m1 * C3 ) * cos ( m1 *t ) + ( -( b ^2 - 3 * m2 ^2) * b * C2 + (3 * b ^2 - m2 ^2) * m2 * C4 ) *cos ( m2 * t) + ( ( m1 ^2 - 3 * b ^2) * m1 * C1 + (3 * m1 ^2 - b ^2) * b* C3 ) * sin ( m1 * t ) + ( -(3 * b ^2 - m2 ^2) * m2 * C2 - ( b ^2 - 3 *m2 ^2) * b * C4 ) * sin ( m2 * t ) ) ;3334sol = [ x1 ( t ) ; x2 ( t ) ; x3 ( t ) ; x4 ( t) ];3536endЛистинг Б.3:findRootsChebFun.m – вспомогательная функция для поисканулей функции, использующая пакет Chebfun.1234function out = findRootsChebFun ( symFunc , tInterval )fun = matlabFunction ( symFunc );out = roots ( chebfun ( fun , tInterval ) ) ;endЛистинг Б.4:determineNextSwitch.m – функция нахождения моментапереключения и соответствующей точки траектории.12function [ tSwitch , xSwitch ] = determineNextSwitch ( odeNum , m1 , m2 ,b , x0 ,tEnd , vpaPrecision )107syms t positivesyms X_1 X_2 X_3 X_4 real345symSolOde = str2func ([ ’ symSolOde ’ , int2str ( odeNum ) ]) ;sol = symSolOde ( m1 , m2 , b ) ;6789101112sol_x1 = subs ( sol (1) ,(3) , x0 (4) } ) ;sol_x2 = subs ( sol (2) ,(3) , x0 (4) } ) ;sol_x3 = subs ( sol (3) ,(3) , x0 (4) } ) ;sol_x4 = subs ( sol (4) ,(3) , x0 (4) } ) ;{ X_1 , X_2 , X_3 , X_4 } , { x0 (1) , x0 (2) , x0{ X_1 , X_2 , X_3 , X_4 } , { x0 (1) , x0 (2) , x0{ X_1 , X_2 , X_3 , X_4 } , { x0 (1) , x0 (2) , x0{ X_1 , X_2 , X_3 , X_4 } , { x0 (1) , x0 (2) , x013% Define the moment of switching numerically :tRoots = findRootsChebFun ( sol_x3 , [0 , double ( tEnd ) ]) ;tRoots = tRoots ( abs ( tRoots ) > 1e -12) ;14151617if ~ isempty ( tRoots )digitsOld = digits ;digits ( vpaPrecision ) ;18192021% Specify the moment of switching using VPA :tSwitch = vpasolve ( sol_x3 , t , tRoots (1) );if tSwitch < 0return ;end222324252627% Define the coordinate of switching using VPA :xSwitch = formula ( ...[ vpa ( subs ( sol_x1 , t , tSwitch ) , vpaPrecision ) ;vpa ( subs ( sol_x2 , t , tSwitch ) , vpaPrecision ) ;vpa ( subs ( sol_x3 , t , tSwitch ) , vpaPrecision ) ;vpa ( subs ( sol_x4 , t , tSwitch ) , vpaPrecision ) ]) ;digits ( digitsOld ) ;28293031323334else35tSwitch = tEnd ; xSwitch = [];36end3738end108Листинг Б.5:integrateTrajectory.m – функция, для моделированиятраектории системы Фиттса с нелинейностью () = sign().1function traj = integrateTrajectory ( m1 , m2 , beta , x0 , tEnd , vpaPresision)234digitsOld = digits ;digits ( vpaPrecision ) ;567a1 = (2*( m1 ^2+ beta ^2) ) * beta +2* beta *( m2 ^2+ beta ^2) ;a0 = ( m1 ^2+ beta ^2) *( m2 ^2+ beta ^2) ;89syms t positive1011traj = { } ; tCurr = 0;1213141516171819202122if x0 (3) < 0[ tSwitch , xSwitch ] = determineNextSwitch (1 , m1 , m2 , beta , x0 ,tEnd - tCurr , vpaPresision ) ;traj = cat (1 , traj , { 1 , x0 , tSwitch } ) ;x0 = xSwitch ;else[ tSwitch , xSwitch ] = determineNextSwitch (2 , m1 , m2 , beta , x0 ,tEnd - tCurr , vpaPresision ) ;traj = cat (1 , traj , { 2 , x0 , tSwitch } ) ;x0 = xSwitch ;endtCurr = tCurr + tSwitch ;2324252627282930313233while tCurr < tEndif x0 (4) > 0[ tSwitch , xSwitch ] = determineNextSwitch (2 , m1 , m2 , beta ,x0 , tEnd - tCurr , vpaPresision ) ;traj = cat (1 , traj , { 2 , x0 , tSwitch } ) ;x0 = xSwitch ;elseif x0 (4) < 0[ tSwitch , xSwitch ] = determineNextSwitch (1 , m1 , m2 , beta ,x0 , tEnd - tCurr , vpaPresision ) ;traj = cat (1 , traj , { 1 , x0 , tSwitch } ) ;x0 = xSwitch ;else1093435363738394041424344454647if a0 * x0 (1) + a1 * x0 (2) > 1[ tSwitch , xSwitch ] = determineNextSwitch (1 , m1 , m2 ,beta , x0 , tEnd - tCurr , vpaPresision ) ;traj = cat (1 , traj , { 1 , x0 , tSwitch }) ;x0 = xSwitch ;elseif a0 * x0 (1) + a1 * x0 (2) < -1[ tSwitch , xSwitch ] = determineNextSwitch (2 , m1 , m2 ,beta , x0 , tEnd - tCurr , vpaPresision ) ;traj = cat (1 , traj , { 2 , x0 , tSwitch }) ;x0 = xSwitch ;elsefprintf ( ’ Sliding mode detected !\ n ’) ;if x0 (2) > 0tSliding = (1 -( a0 * x0 (1) + a1 * x0 (2) ) ) / a0 / x0 (2) ;traj = cat (1 , traj , { 3 , x0 , tSliding } ) ;x0 = [(1 - a1 * x0 (2) ) / a0 ; x0 (2) ; 0; 0];48traj = cat (1 , traj , { 2 , x0 , tSwitch } ) ;x0 = xSwitch ;elseif x0 (2) < 0tSliding = ( -1 -( a0 * x0 (1) + a1 * x0 (2) ) ) / a0 / x0 (2) ;traj = cat (1 , traj , { 3 , x0 , tSliding } ) ;x0 = [( -1 - a1 * x0 (2) ) / a0 ; x0 (2) ; 0; 0];49505152535455traj = cat (2 , traj , { 2 , x0 , tSwitch } ) ;x0 = xSwitch ;5657else58fprintf ( ’ Equilibrium interval is reached !\ n ’);break ;5960end61end62endtCurr = tCurr + tSwitch ;63646566endsave ( ’ traj .