Решение обыкновенных дифференциальных уравнений в среде Scilab

Солвер служит для решения обыкновенного дифференциального уравнения (ОДУ)

dy/dt=f(t,y) , y(t0)=y0

Синтаксис

y=ode(y0,t0,t,f)

[y,w,iw]=ode([type],y0,t0,t [,rtol [,atol]],f [,jac] [,w,iw])

[y,rd,w,iw]=ode("root",y0,t0,t [,rtol [,atol]],f [,jac],ng,g [,w,iw])

y=ode("discrete",y0,k0,kvect,f)

Параметры

y0 : действительное число или матрица (условия инициализации)

t0 : действительный скаляр (время инициализации)

t : действительный вектор. Задает интервал времени для которого вычисляется решение уравнения.

f : внешний параметр (функция или строка или список).

type : строковая переменная принимающая одну из следующих значений: "adams" "stiff" "rk" "rkf" "fix" "discrete" "roots". adams применяют при решении дифференциальных уравнений или систем методом прогноза-коррекции Адамса; stiff указывают при решении жестких задач; rk используют при решении дифференциальных уравнений или систем методом Рунге-Кутта четвертого порядка; rkf указывают при выборе пятиэтапного метода Рунге-Кутта четвертого порядка; fix тот же метод Рунге-Кутта, но с фиксированным шагом;

rtol, atol : действительные константы и действительные константы того же размера, что и y, по умолчанию rtol=0.00001, atol=0.0000001, при использовании параметров rkf и fix - rtol=0.001, atol=0.0001;

jac : матрица, представляющая собой якобиан правой части жесткой системы дифференциальных уравнений, задают матрицу в виде внешней функции вида J=jak(t,y);

w, iw : векторы, предназначенные для сохранения информации о параметрах интегрирования, которые применяют для того, чтобы последующие вычисления выполнялись с теми же параметрами.

ng : целое число.

g : внешний параметр (функция или строка или список).

k0 : целое число (начальное время).

kvect : целочисленный вектор.

Листинг

function yd=f(t,x),yd=-x+sin(t*x),endfunction;

x0=1.5;

t0=0;

t=0:1:35;

y=ode(x0,t0,t,f);

plot(t,y)

Листинг

//Функция, описывающая систему дифференциальных уравнений

function dy=syst(t,y)

dy=zeros(2,1);

dy(1)=cos(y(1)*y(2));

dy(2)=sin(y(1)+y(2)*t);

endfunction

//Решение системы дифференциальных уравнений

x0=[0;0];

t0=0;

t=0:0.1:10;

y=ode(x0,t0,t,syst);

//Формирование графического решения

plot(t,y)

xtitle('diff equation','t, c','y')

legend('y')

Листинг

function dx=syst2(t,x) //Функция задающая систему ОДУ

dx=zeros(3,1);

dx(1)=-7*x(1)+7*x(2);

dx(2)=157*x(1)+x(2)-1.15*x(1)*x(3);

dx(3)=0.96*x(1)*x(2)-8.36*x(3);

endfunction

//Решение ОДУ

x0=[-1;0;1];

t0=0;

t=0:0.01:2;

y=ode("stiff",x0,t0,t,syst2);

plot(t,y);

xgrid();

Листинг

function F=FF(t,x)

F=[-4*x(1)-13*x(2)+exp(t); x(1)];

endfunction

//Решение системы дифференциальных уравнений

X0=[1;-1];

t0=0.25;

t=0.25:0.05:2;

y=ode("stiff",X0,t0,t,FF);

//Вывод графика решения

plot(t,y);

xgrid();