Mоделирование процессов и систем в Matlab (966709), страница 21
Текст из файла (страница 21)
Ф мойя - строка с иненен процедуры. к которой $ обращается процедура тмп; $ С вЂ” текущий нанент вренени: д у — вектор текуцнх значений перененных состояния; т т — вычисленные значения производных т(1) - Оу(т)/ОС: т Л вЂ” шаг интегрирования; т С вЂ” предыдущий нонент вренени; д у — предыдущее значение вектора перененных состояния; д Выходные паранетры: д Самс — новый нонент вренени: д уомт — новое значение вектора перененных состояния через шаг интегрирования $ Расчет проненуточных значений производных К1 - теча)(2ртмп.мртмп,т.у); )г2 " (еча)(2рйю.мртип С + М/З.у + П/Зяк1): КЗ - теча)(2ртып,мртип,т и 2ЯП/З,у + ЛЯД2 - П/ЗЯД1): К4 - теча)(2ртмп.мрйаЛ+П.у + Пя(КЗ + К1 - К2)): $ Расчет новых значений вектора перененных состояния Спит - С + М; усттт у + Пя(К) + 3'ЧС2 + ЗяКЗ + К4)/В: $ Конец процедуры ЙК043ш Такая форма представления процедуры вычисления правых частей дифференциальных уравнений неудобна.
Во-цервых, процедуру вида ЕМ1 нельзя использовать при,интегрировании процедурами МАТЮКАВ е1е23 ы ООВ45 (последние требуют, Урок 2 ° Програииирование в среде МАТО(В чтобы в процедуре вычисления правых частей уравнения было только два вхолных параметра, а в процедуре ГМ1 их три).
Во-вторых, такая форма вызовет необходимость создания новых М-файлов методов численного интегрирования. Этого можно избежать, представив имя дополнительной функции Мруцп как глобальную переменную. Тогда процедура вычисления правых частей уравнения может быть записана так: Рцпсттоп г - ГМ2(г.у): $ Процедура правых частей уравнения физического наятника Г Осуществляет расчет вектора з" $ произвадних вектора "у" перененних состояния по форнулан: з гП) - у(2): т з(2) = -зтп(у(1)) + В(г,у) 3 Входные паранетри: $ щр(цп — иня процедури з((,у) — глобальная перененная $ щргцп - 'ь'; $ 1 — текущий нонент вренени; Х у — текущее значение вектора перененних состояния: $ Виходние паранетри: д з — вектор значений производних от перененних состояния Я1ооа! КРГР г(1) у(2)„.
г(2) - -Юп(у(1)) + (еча1(МРГОМЛ.у): $ Конец процедуры ГМ2 Теперь процедура ГМ2 имеет только два входных параметра, передаваемых через заголовок, и может быть использована любой процелурой численного лтетода интегрирования, включая процедуры от(е23 и о()е45. Необходимо лишь помнить, что в основной программе переменной МРГ()М надо присвоить некоторое символьное значение (нмя функции, которая будет использована в процедуре вычисления правых частей уравнения), и эта переменная должна быть объявлена как глобальная.
Например, если будет использована Ранее созданная процедура МощГцп1, в файле-сценарии должны присутствовать такие строки: р)оца) МРГВМ КРГР 'Мотрщ1 ': Программа моделирования движения маятнина Ранее были рассмотрены проблемы, с которыми приходится сталкиваться во время создания сложных программ, и средства их преодоления. Теперь на основе полученных знаний попробуем составить и испытать в работе одну из довольно сложных комплексных программ.
Создадим программу, которая позволила бы моделировать движение физического маятника с вибрирующей точкой полвеса путем численного интегрирования дифференциального уравнения этого движения. Дифференциальное уравнение движения маятника для этой задачи можно принять таким: .7 тр+ Я(р+ лгЩ1+ и зш(оз(+ е„))зшгр = -тля(п зш(о)с+ в )соз тр, 11з Программа моделирования движения маятника где.т' — момент инерции маятника; гг — коэффициент демпфирования; гн4 — коэффициент жесткости; 脄— амплитуда виброперегрузки точки подвеса маятника в вертикальном направлении; и — амплитуда виброперегрузки в горизонтальном направлении; «Р — угол отклонения маятника от вертикали; от — частота колебаний точки подвеса; е е„— начальные фазы колебаний точки подвеса в горизонтальном и вертикальном направлениях.
Нужно создать такую программу, которая позволяла бы вычислять закон изменения угла отклонения маятника от вертикали во времени при произвольных, устанавливаемых пользователем значениях всех вышеуказанных параметров маятника и поступательного движения основания, а также при произвольных начальных условиях. Вычисления будем осуществлять путем численного интегрирования с помощью стандартной процедуры обе45. Преобразование уравнения Для подготовки дифференциальных уравнений к численному интегрированию прежде всего необходимо привести эти уравнения к нормальной форме Коши. Желательно также представить их в безразмерной форме. Для представленного уравнения это было сделано в предыдущем разделе. И-файл процедуры вычисления правых частей уравнений Следующим шагом подготовки программы является написание и запись на диск текста процедуры вычисления правых частей полученной системы ОДУ в форме Коши.
Эта процедура была создана в предыдущем разделе и в окончательном варианте она имеет такой вид: Файл ГМ2.ы' Гипс(1 оп г ЕИ2(т.. т): $ Процедура правих частей уравнения Физического наятника $ Осуществляет расчет вектора "г" производных вектора д "у" перененных состояния по Форнулан: $ г(1) - у(2): $ г(2) " -з1п(уП)) + 5(т.у) $ Входные параметры: 3 С вЂ” текущий момент вренени: $ у — текущее значение вектора перененних состояния: $ ИРРОИ вЂ” имя процедуры 5(т.у) — глобальная переменная д ИРГОИ - '5': $ Выходные параметри: $ г — вектор значений производных от перененных состояния В1оЬа1 ИРГОЙ г(1) - у(2): г(2) - -з1п(у(1)) + тета)(ИРРОИЛ .у): а 7„' т Конец процедуры РИ2 Исходныс тексты всех описывясмых в книге нрогрзмм можно загрузить с ввб-сервере издательства, Урок 2 ° Программирование в среде МАТСАВ ~ ПРИМЕЧАНИЕ Обратите внимание на неболывое, но существенное отличие данной процедуры от аналогичной процедуры, приведенной в предыдущем разделе, — наличие в конце процедуры операции транспонирования вектора производных.
Это обусловлено тен, что процедура аде45 требует, чтобы вектор производных обязательно был представлен в виде столбца. В качестве дополнительной процедуры, используемой в процедуре ЕМ2, выберем ранее созданную процедуру МсвГМ1, которую запишем в файл Моти ГМ1.
Файл Мое ГМ1.в Уопсс)оп в - МояГМ1П.у): $ Вычисление Монентов Сил, действуюкнх на физический наятник Д Осуиествляет расчет нонента "в" снл Г по форнуле: д в --2ябгяу(2) - (тахьз1п(пьяС + ех)"соя(уП)) +... д + гду йтп(пи"С + еу)*зтп(уП)). Х Козй)нцненты пцзедаются в процедуру через $ глобальный вектор КМ1 - (бг.пюу.вяхлм.еу.ех] В1о)в) КМ1 в = -2яКМ1П)'у(2) - КМИЗ)*зтп(КМ1(4)'"С + КМН6))*соя(уП)) -.. КМ1(2) з1п(КМ1(4)"С + КМ1(5))"Ып(уП)): Д Конец процедуры МогГМ1 Очевидно, что в вызывающем файле-сценарии надо предусмотреть объявление имени дополнительного файла МогпЕМ1 как глобальной переменной МРГОМ, а также обеспечить объявление глобальной переменной по имени КМ1 и определение значений зтого числового массива из пяти элементов.
Управляющий файл-сценарий Главный файл создадим в соответствии с рекомендациями, привеленными в предыдущем разделе. Файл Язмауатп2.щ д Г1П(ауатл2 д упрааляююая програнна исследования движения физического наятннка Д установЛеннОГО на ПоСтупатЕЛьна Внбрируояен основании Г1знауасп2 ?азсауйа К - вело('Что делать? ', 'Продолжить работу'. 'Закончить работу' ): гб К нлт)е К 1 Ягнауакп2 Мега Е1гтцуасп2 Уабго К - веты('Что делать? '. 'Продопяить работу'.
'Закончить работу'): епб вх) с1еаг В1оба1 с1еаг $ Конец Е1гИауасп2 Как видим, программа лишь вызывает три дополнительных (райла-сценария— ЯгМауаП)2 2аякачйа, Я2МауаСП2 Мепц и Я2МауаСП2 уадго. 1]оэтому нужно еще создать зти М-файлы. Программа моделирования движения маятника 115 Файл-сценарий заставки Мы уже отмечали, что файл-сценарий заставки должен содержать операторы вывода на экран информации об основных особенностях математической модели, реализованной в программе, и операторы ввода исходных значений параметров этой модели. Ниже приведен текст М-файла Вхйауатп2 2аз(аы)га.
Файл Мяиауасп2 2азтавйа.гл $ Г1зизуавп2 Хазтатда К Часть (вывод заставки на »кран) програнны Г!гпауатп2 д Ввод "вшитых" значений ЗРгоргаш - 'Е!зиаутп2.н'; зпмге -'Лазарев Ю.ф. '. КМ1 - (О 0 0 0 0 0); КРЕОН " 'МошГш1': 0)оба) КИ1 ИРЕОИ СГ1па) - 2"р!"5: Г10 " р1(180: Г!СО " 0; с)с одзрн 'Это програнна, осуществляовая интегрирование уравнения':... 'Визическаго Иаятника при поступательной вибрации точки подвеса". 'в форне':... ' Г!" + з1п(Г1) " -2"ба*К!'' - ';... ' - мяу»в!п(пы»С + еу)яз!п(Г!) - пахяз1п[пц»С + ех)»соз((!)':.
'где Г! — угол отклонения наятника от вертикали, '.-... ' бг — относительный козффициент затухания. ':... ' пц — относительная частота вибрации точки падвеса, ';... ' пгзг. пвх — анплитуды виброперегрузкн в вертикальнон ';... ' н гаризонтальнон направлениях соответственно. ';... ' еу.ех — начальные фазы колебаний в вертикальнон ';... ' и горизонтальнан направлениях соответственно„ '.-.. ' КМ1 - (ба,пшу,пшх.пц.еу.ех) — натрица козффициентав')) К Конец Г!змауатп2 Хазтавйа В этом файле осуществляется присваивание исходных (лвшитых») значений всем параметрам заданного дифференциального уравнения, а также параметрам численного интегрирования — начальным условиям движения маятника и длительности процесса интегрирования.
Часть этих параметров объединяется в единый глобальный вектор КИ1. Одновременно переменной МРГ()И, которая будет использоваться при интегрировании, присваивается значение 'ИопЕи1 '. Файл меню Содержимое файла меню гтдМауа(п2 Мепц приведено ниже. Файл йяМауатл2 Мепц.гп К Е1зпауатп2 Мепц К Чань (осуществлявшая диалоговое изненение данных) К програины Е1гмауасп2 К 1: ьн!)е К < 10 0150(' ') б!зр('Сейчас установлено') 01зр((зргтпСГ('Начальный угол (градусы) - ДО', Г10"100/01)... 11б Урон 2 ° Программирование в среде МАТСАВ зрг1пСП 'Начальная скорость - $9'. Г!СОН) 01зр<зрг!пС((' Число периодов = гд'.