Maple 9.5/10 в математике, физике и образовании

11.1.8. Анализ и моделирование линейных систем операторным методом

Произвольные линейные системы могут анализироваться и моделироваться хорошо известным (особенно в электротехнике и радиотехнике) операторным методом. При этом методе система и ее воздействие представляются операторными выражениями, т. е. в виде функций параметра — оператора Лапласа s (в литературе встречается и обозначение p). Не вникая в детали этого общеизвестного метода, рассмотрим конкретный пример (файл linsys). Он, для сравнения с предшествующими примерами, дан для системы второго порядка, хотя в данном случае никаких ограничений на порядок системы нет.

Для начала зададим инициализацию применяемых пакетов расширения

> restart:with(plots): readlib(spline): with(inttrans):
Warning, the name changecoords has been redefined

Далее зададим операторные выражения для коэффициента передачи системы G и входного сигнала R (в виде единичного перепада) и вычислим с упрощением их произведение:

> G := K/(M*s^2+C*s+1); R := 1/s;

> X := simplify(R*G);

Теперь, используя обратное преобразование Лапласа, найдем временную зависимость реакции системы в аналитическом (что наиболее ценно) виде:

> h := simplify(invlaplace(X, s, t));

Теперь мы можем построить график этой зависимости для конкретных значений М, С и K:

> h1 : = subs(M=1,C=0.75,K=1,h);

h1 := 0.5393598900 I (-1.854019622 I + 1.854049622 I e^{(-0.3750000000 t)} cosh (0.9270248110 I t) + 0.75 e^{(-0.3750000000 t)} sinh(0.9270248110 I t))

> linresp := plot(h1, t=0..20, axes=boxed, color=black): display(linresp);

Вид этой зависимости представлен на рис. 11.19. Он соответствует реакции системы второго порядка для случая затухающих колебаний.

Рис. 11.19. Одна из временных зависимостей реакции системы второго порядка

А теперь зададимся целью наглядно проиллюстрировать изменение временной зависимости реакции системы при изменении параметра С от 0 до 2 при М=1 и K=1. Для этого выполним следующие вполне очевидные команды:

> х := subs(M=1, K=1, h);

> plot3d(x, С=0..2, t=0..20, axes=boxed);

Соответствующий график показан на рис. 11.20. Он прекрасно иллюстрирует переход от апериодического режима при С=2 к колебательному при С= 0 при изменении времени от 0 до 20.

Рис. 11.20. Динамика развития колебаний в системе при изменении параметра С

Аналогичным образом можно построить трехмерный образ временной зависимости реакции системы для М=1, С=0.25 и изменении параметра K от 0 до 3. Для этого надо выполнить команды:

> x1 := subs(М=1, С=0.25, xt);

x1 := 0.5039526307 I K (-1.984313483 I + 1.984313483 I e^{(-0.1250000000 t)} cosh (0.9921567415 I t) + 0.25 е^{(-0.1250000000 t)} sinh(0.9921567415 I t))

> plot3d(x1, K=0..3, t=0..20, axes=boxed);

Диаграмма временных зависимостей представлена на рис. 11.21.

Рис. 11.21. Динамика развития колебаний в системе при изменении параметра K

Представленные на рис. 11.20 и 11.21 диаграммы дают весьма наглядное представление о динамике поведения рассмотренной системы. Но еще важнее то, что просто изменением операторной записи G и R по описанной методике можно анализировать и наглядно представлять работу множества линейных систем.

11.2. Моделирование динамических задач и систем

11.2.1. Расчет траектории камня с учетом сопротивления воздуха

Вы хотите метнуть камень в огород вашего вредного соседа? Разумеется во время его отсутствия. Давайте промоделируем эту ситуацию, предположив два актуальных случая: дело происходит на Земле в условиях, когда наша планета лишилась воздуха и когда, слава богу, он все же есть. В первом случае сопротивления воздуха нет, а в другом сопротивление воздуха есть и его надо учитывать. Иначе камень упадет в ваш огород, а не в огород соседа!

Учет сопротивления воздуха не просто усложняет задачу нашу задачу. Он делает ее нелинейной. В связи с этим мы применим численные методы решения дифференциальных уравнений. Кроме того, учитывая громоздкость документов, описывающих приведенные ниже задачи, перейдем к их записи прямо в тексте книги.

Итак, пусть подвернувшиеся под руку камни с массой 500 и 100 грамм брошены под углом 45 к горизонту со скоростью V_O=20 м/с. Найдем их баллистические траектории, если сила сопротивления воздуха Fmp=A*V, где А=0,1 Н∙с/м. Сравним их с траекториями, получающейся без учета сопротивления воздуха. Документ с решением этой задачи, описанным ниже, представлен в файле balist.

Начнем с подключения пакета plots, нужного для визуализации данной задачи:

> restart; with(plots):
Warning, the name changecoords has been redefined

Составим параметрические уравнения для проекций скорости на оси координат:

> Vox:=Vc*cos(alpha);Voy:=Vo*sin(alpha);

Vox:= Vo cos(α)

Voy:= Vo sin(α)

Мы рассматриваем два случая: камень массой 500 г и камень массой 100 г. Поскольку для каждого случая мы предусматриваем расчет в двух вариантах (с учетом сопротивления воздуха и без такого учета), то мы должны составить 4 системы дифференциальных уравнений (ДУ). Каждая система состоит из двух ДУ второго порядка и вид этих систем известен из курса физики. Ниже представлено задание этих систем ДУ (для первой системы дан вывод ее вида):

> sys1:=massa[1]-diff(x(t), t$2) =
-A[1]*diff(x(t),t),massa[1]*diff(y(t), t$2) = -A[1]*(diff(y(t),t))-massa[1]*g;