Решения систем ОДУ

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

Уравнение нужно переписать так, чтобы самая большая производная осталась одна в левой части:

Теперь определим массу тела, жесткость пружины, коэффициент затухания и вынуждающую силу:

m = 2.5

C = 1.4e3

b = 2.3

A = 5.2

omega = 50*pi

f(x) = A*sin(omega*x)

Теперь начните создавать блок ОДУ с введением команды ode(имя, переменная), где имя - это идентификатор блока ОДУ (задайте любое имя), и переменная - это переменная интегрирования. Выполните:

ode(s1, t)

далее введите уравнение:

x'' = (1/m)*(-C*x - b*x' + f(x))

как видите, для обозначения производной используется апостроф.

потом введите команду end, чтобы закончить создание блока ОДУ:

end

создайте вектор, который содержит значения переменной интегрирования:

t = < 0: .001 : 5 >

Создайте вектор с начальными условиями. Поскольку в блоке одно дифф. уравнение 2-ого порядка, вектор с начальными условиями должен содержать два значения (начальное перемещение и начальная скорость). Выполните:

X0 = .12

V0 = -1.5

x0 = < X0, V0 >

и, наконец, используйте команду odesolve(name, var, x0) для решения блока ОДУ:

R = odesolve(s1, t, x0)

Теперь R содержит все посчитанные значения искомых функций в таблице (R - это 2D массив). Первая колонна - это вектор переменной интегрировария t, вторая колонна содержит перемещения x, и третья колонна содержит скорости v. Теперь можно извлечь, скажем, 2-ую колонну следующим образом:

x = vector(R[ ,2])

С помощью команды explode(имя ОДУ, имя Массива), можно автоматически создать все необходимые векторы с данными, чтобы для каждого вектора не вводить команду vector(R[ ,n])

explode(s1, R)

и, наконец, отобразите результаты:

plot(t, x)

В результате рабочий лист должен выглядеть примерно так:

Если нужно отредактировать уравнение в блоке, введите команду modify(имя ОДУ), затем введите нужное уравнение заново и выйдите из блока с помощью команды end.