Учебное пособие для студентов высших учебных заведений
Скачать 5.41 Mb.
|
function z = FM2(t,y); % Процедура правых частей уравнения Физического Маятника. % Осуществляет расчет вектора “z” % производных вектору "y" переменных состояния по формулам: % z(1)=y(2); % z(2)= - sin(y(1)) +S(t,y), % Входные параметры: % mpfun - имя процедуры S(t,y) - глобальная переменная % mpfun = 'S'; % t - текущий момент времени; % y - текущее значение вектора переменных состояния; 2.6. Создание функций от функций 124 % Выходные параметры: % z - вектор значений производных от переменных состояния. global MPFUN z(1) = y(2); z(2) = - sin(y(1)) + feval(MPFUN,t,y); % Конец процедуры FM2 Теперь процедура FM2 имеет только два входных параметра, передаваемых через заголовок, и может быть использована любой процедурой численного мето- да интегрирования, в том числе - процедурами ode23 и ode45. Необходимо лишь помнить, что в основной программе переменной MPFUN надо присвоить некото- рое символьное значение (имя функции, которая будет использована в процедуре правых частей), и она должна быть объявлена как глобальная. Например, если бу- дет использована ранее созданная процедура MomFun1, в Script-файле должны присутствовать строки global MPFUN MPFUN = ‘MomFm1'; 2.6.3. Задания Задания 2.1 - 2.13. Создайте М-файл метода численного интегрирования дифференциальных уравнений в соответствии с формулами, приведенными в таб- лицах 2.1 и 2.2. Таблица 2.1. Методы Рунге-Кутта y m+1 = y m + h F(t m ; y m ) N% вар Формула метода Вспомогательные величи- ны Название методу 1 F=k1 k1=Z(t m ;y m ) Ейлера 2 F=(k1+k2)/2 k1=Z(t m ;y m ); k2=Z(t m +h;y m +hk1) модифициро- ванный Ейлера 3 F=Z(t m +h/2; y m +hk1/2) k1=Z(t m ;y m ) 4 F=(k1+4k2+k3)/6 k1=Z(t m ;y m ); k2=Z(t m +h/2;y m +hk1/2); k3=Z(t m +h;y m +h(2k2-k1)) Хойне 5 F=(k1+3k3)/4 k1=Z(t m ;y m ); k2=Z(t m +h/3;y m +hk1/3); k3=Z(t m +2h/3;y m +2hk2/3) 6 F=(k1+2k2+2k3+ +k4)/6 k1=Z(t m ;y m ); k2=Z(t m +h/2;y m +hk1/2); k3=Z(t m +h/2;y m +hk2/2); k4=Z(t m +h;y m +hk3) Рунге-Кутта 2.6. Создание функций от функций 125 7 F=(k1+3k2+3k3+ +k4)/8 k1=Z(t m ;y m ); k2=Z(t m +h/3;y m +hk1/3); k3=Z(t m +2h/3;y m +h(k2- k1/3)); k4=Z(t m +h;y m +h(k1-- k2+k3)) Таблица 2.2. Многошаговые методы N% вар Формула прогнозу Формула коррекции Название методу 8 y m+1 =y m-1 +2h(t m ;y m ) y m+1 =y m +h[Z(t m+1 ;y * m+1 )+ +Z(t m ;y m )]/2 9 y m+1 =y m +h[Z(t m ;y m )- -Z(t m-1 ;y m-1 )]/2 y m+1 =y m +h[Z(t m+1 ;y * m+1 )+ +Z(t m ;y m )]/2 10 y m+1 =y m +h[23Z(t m ;y m )- -16Z(t m-1 ;y m-1 )+ +5Z(t m-2 ;y m-2 ]/12 y m+1 =y m +h[5Z(t m+1 ;y * m+1 )+ +8Z(t m ;y m )- -Z(t m-1 ;y m-1 )]/12 11 y m+1 =y m +h[55Z(t m ;y m )- -59Z(t m-1 ;y m-1 )+ +37Z(t m-2 ;y m-2 )- -9Z(t m-3 ;y m-3 )]/24 y m+1 =y m +h[9Z(t m+1 ;y * m+1 )+ +19Z(t m ;y m )- -5Z(t m-1 ;y m-1 )+ +Z(t m-2 ;y m-2 )]/24 Адамса- Башфорта 12 y m+1 =y m-3 +4h[2Z(t m ;y m )- -Z(t m-1 ;y m-1 )+ +2Z(t m-2 ;y m-2 )]/3 y m+1 =y m-1 +h[Z(t m+1 ;y * m+1 )+ +4Z(t m ;y m )+ +Z(t m-1 ;y m-1 )]/3 Мілна 13 y m+1 =y m-3 + +4h[2Z(t m ;y m )- -Z(t m-1 ;y m-1 )+ +2Z(t m-2 ;y m-2 )]/3 y m+1 ={9y m - y m-2 + +3h[Z(t m+1 ;y * m+1 )+ +2Z(t m ;y m )- -Z(t m-1 ;y m-1 )]}/8 Хеммінга Задание 2.14. Создайте М-файл процедуры правых частей дифференциаль- ных уравнений движения двухстепенного гироскопического компаса: , 0 ) sin ) ( cos ) ( sin cos ( * * )] sin ) ( cos ) ( cos cos ( [ 3 3 2 1 = + + − + − + β β β ϕ ω β β β ϕ ω β t u t u t u t u J H J N E E N && где - моменты инерции гирокомпаса; J J 1 2 , H - его собственный кинетический момент; β - угол отклонения главной оси гирокомпаса от плоскости географиче- ского меридиана места ; ϕ - географическая широта места объекта, на котором установлен гирокомпас; ), sin( ) ( N Nm N t u t u ε ω + = ) sin( ) ( E Em E t u t u ε ω + = - соответственно северная и восточная составляющие угловой скорости поворота основания, на котором установлен гирокомпас. 2.6. Создание функций от функций 126 Задание 2.15. Создайте процедуру правых частей дифференциальных урав- нений, которые описывают динамику объемов популяций хищников и жертв и известны как модель Вольтерра: x t 1 ( ) x t 2 ( ) ; 2 1 21 1 22 2 2 1 12 1 11 1 x x a x a x x x a x a x − = + − = & & Задание 2.16. Создайте процедуру правых частей дифференциального уравнения углового движения торпеды в горизонтальной плоскости, которая управляется нелинейным исполнительным элементом: , 0 ) ( 2 2 = + + ψ ψ ψ kF dt d R dt d J Процедура должна предусматривать возможность использования нескольких су- щественно нелинейных законов управления ) (x F , в частности, релейного с зона- ми нечувстввительности и гистерезисом: если , то F( &x > 0 x c x b b x b c x b ) , , ; = − < − − < < + > ⎧ ⎨ ⎪ ⎩⎪ при при при 1 1 2 2 0 если же , то &x < 0 F x c x b b x b c x b ( ) , , = + > − < < − < ⎧ ⎨ ⎪ ⎩⎪ при при при 1 2 1 2 0 Задание 2.17. Создайте процедуру правых частей дифференциального уравнения углового движения искусственного спутника Земли, управляемого по логическому закону: ); , ( ϕ ω ω Φ = k dt d J ω ϕ = dt d , где ω - угловая скорость движения спутника; J - его момент инерции; ) , ( ϕ ω Φ - заданная логическая нелинейная функция, которую определим с помощью такой таблицы: Значение функции ) , ( ϕ ω Φ Знак ϕ Знак ω - 0 + - +1 0 0 0 +1 0 -1 + 0 0 -1 Задание 2.18. Создайте процедуру правых частей дифференциальных уравнений движения волчка со сферическим подпятником, установленным в ко- нической лунке: 2.6. Создание функций от функций 127 ⎪ ⎪ ⎪ ⎪ ⎩ ⎪ ⎪ ⎪ ⎪ ⎨ ⎧ + = + − = + = + = = , cos sin , sin sin sin cos cos , sin , cos sin , 1 1 2 2 1 2 1 1 1 р 2 р 2 1 р δ δ δ δ δ δ δ δ δ δ δ δ ς ξ ς ξ η ς ς η η ξ ξ K K J K K K J M mgl dt dK M mgl dt dK M dt dK e e & & где - экваториальный момент инерции волчка относительно точки опоры; J e H - кинетический момент волчка; δ δ 1 2 , - углы отклонения оси волчка от вертикали; mgl - опорный маятниковый момент волчка; ς - составляю- щиеьмомента сил трения в подпятнике волчка; - проекции кинетиче- ского момента волка на неподвижные оси. Моменты сил трения можно представить следующими зависимостями: M M M т т т р р , , ξ η р ς р ς K K K ξ η , , M M M т т т р р , , ξ η ], 2 / ) (sin 1 [ sin cos }; 2 / ) (sin sin cos ] 2 / ) (sin 1 [ {sin ; cos cos cos 2 1 2 р 2 1 2 2 2 р 2 1 2 р α δ δ α δ δ α δ δ δ α ς η ξ − = + − − = − = C M C M C M где ), ( cos 1 H sign B mgl l R k C α = а 2 / ) cos sin (sin sin cos cos cos 1 2 2 1 2 2 2 2 2 2 1 2 2 δ δ δ α δ δ α + − − = B Выше использованы обозначения: k - коэффициент трения материала под- пятника и материала опоры; R - радиус сферы подпятника; α - угол между обра- зующей конуса лунки и плоскостью горизонта. Взаимосвязь проекций с собственным кинетическим моментом K K K ξ η , , ς H и другими кинематическими величинами определяется соотношениями: ; sin cos cos sin cos cos 1 2 1 2 2 1 1 2 δ δ δ δ δ δ δ δ ξ & & e e J J H K + − = ; cos sin 2 2 1 2 δ δ δ η & e J H K + = cos sin cos sin sin cos 1 2 1 2 2 1 1 2 δ δ δ δ δ δ δ δ ς & & e e J J H K + + − = Задание 2.19. Создайте процедуру правых частей дифференциальных уравнений гироскопа в кардановом подвесе (ГКП), установленного на неподвиж- ном основании: 2.6. Создание функций от функций 128 ⎪ ⎪ ⎪ ⎩ ⎪ ⎪ ⎪ ⎨ ⎧ + + = + + + + − = − + + + − + + + − = = + − + ), sin( ), sin( cos cos sin , sin )] sin( [ ) sin( cos cos sin 2 ) cos ( 0 0 2 2 3 0 0 2 2 2 2 1 R m L m R m N m t R R dt dH t L L f H J J t R R t N N f H J J J ε ω ε ω β β α β β α β β ε ω ε ω α β β β β β α α β & & & && & & & & && где ; 1 2 1 Z X J J J + = ; 1 1 2 Z e X J J J J − + = ; 1 3 e Y J J J + = - момент инерции внешней рамки карданового подвеса относительно наружной оси подвеса; - моменты инерции внутренней рамки относительно указанных осей; - экваториальный момент инерции ротора гироскопа; X J 2 Z Y X J J J 1 1 1 , , e J β α , - углы поворота главной оси ГКП вокруг наружной и внутренней осей подвеса; H - собственный кинетический момент ГКП; - коэффициенты вязкого трения по внутренней и наружной осям подвеса; - постоянные составляющие моментов внешних сил, направленных по наружной, внутренней осям подвеса и главной оси гироскопа соответственно; - амплитуды гармонических составляю- щих моментов сил, действующих по соответствующим осям; 2 1 , f f 0 0 0 , , R L N m m m R L N , , ω - частота измене- ния гармонических составляющих моментов сил; R L N ε ε ε , , - начальные фазы гармонических составляющих моментов сил. Задание 2.20. Создайте процедуру правых частей дифференциального уравнения гироскопического тахометра (ГТ), установленного на вращающейся основе: , / ] ) / ( ) / ( [ 0 1 1 2 1 1 1 1 H c M u u H J u H J u c x c x f x J Z X Y Z + + − − = + + & & && где - выходной сигнал ГТ; x H - собственный кинетический момент ГТ; - моменты инерции ГТ; - угловая жесткость упругой связи ГТ с основанием; 2 1 , J J c f - коэффициент углового демпфирования; - проекции угловой скорости основания на оси, связанные с ГТ. Последние связаны с проекциями на оси, свя- занные с основанием, соотношениями: 1 1 1 , , Z Y X u u u , ; sin cos ; sin cos 1 1 1 Y Y X Z Z Z X X u u u u u u u u = + = − = β β β β где x c H − = β Проекции угловой скорости основания на оси, связанные с тем же основа- нием, полагать изменяющимися со времени по законам: ). sin( ); sin( ); sin( 0 0 0 Z Zm Z Z Y Ym Y Y X Xm X X t u u u t u u u t u u u ε ω ε ω ε ω + + = + + = + + = 2.6. Создание функций от функций 129 Задание 2.21. Следящая система состоит из задающего элемента, который задает 1 ϑ угол, на который должен повернуться выходной вал следящей системы (вал электродвигателя), формирующего элемента (сельсина), который сравнивает этот угол с углом поворота 2 ϑ выходного вала электрического двигателя и фор- мирует электрический сигнал, пропорциональный синусу разности этих двух уг- лов: ). sin( 2 1 1 1 ϑ ϑ − = m U u Этот сигнал суммируется с сигналом тахогенератора на валу двигателя: ; 1 2 Ђ Ђ k k k u u u u ω = − = Сигнал подается на усилительное устройство, которое представляет со- бой трехпозиционное реле с гистерезисом. Последнее формирует напряжение u в соответствия с зависимостью: u 2 Д u f u z sign u п и x п и u x Д b b a = = > < ⎧ ⎨ ⎩ ( ) ( ) р | ; р | 2 2 2 0 |u | 2 Вращательное движение вала двигателя описывается дифференциальными уравнениями: ⎪ ⎩ ⎪ ⎨ ⎧ = = + ; 2 Ђ Ђ Ђ Ђ Ђ Ђ dt d u k dt d T ω ϑ ω ω Создайте в форме М-файла процедуру вычисления правых частей диффе- ренциальных уравнений следящей системы, считая выходными величинами угол 2 1 ϑ ϑ ϑ − = рассогласования и скорость его изменения. Задание 2.22. Составьте процедуру отыскания точных решений системы линейных однородных дифференциальных уравнений по заданной матрице системы ДУ в форме Коши: A d dt y A y = ⋅ и заданному вектору y0 начальных условий. Задание 2.23. Создайте процедуру правых частей дифференциальных уравнений движения в пространстве трех гравитирующих материальных точек (задача трех тел в небесной механике) ⎪ ⎪ ⎪ ⎪ ⎩ ⎪⎪ ⎪ ⎪ ⎨ ⎧ + − = − − = − − = ), ( 1 ), ( ), ( 3 2 1 3 13 1 3 32 2 2 2 3 32 3 3 21 1 2 2 3 2 1 13 32 3 32 21 2 R R R R R R R R R m m m R m R m dt d R m R m dt d γ γ |