Page 130 - 363_

131

варіант FM1:
function z = FMl(mpfun, t, y);
% Процедура правих частин рівняння Фізичного Маятника.
% Здійснює розрахунок вектора “z” похідних вектора “у” змінних стану
% за формулами:
% z(l)=y(2);
% z(2)=-sin(y(l))+S(t,y),
% де вхідні параметри:
% mpfun - ім'я процедури S(t,y)
% mpfun = 'S';
% t - поточний момент часу;
% у - поточне значення вектору змінних стану;
% вихідні параметри:
% z - вектор значень похідних від змінних стану.
z(l) = y(2);
z(2) = - sin(y(l)) + feval(mpfun,t,y);
% Кінець процедури FM1
Через те, що вид звернення до процедури правих частин змінився,
необхідно також перебудувати і процедуру чисельного методу. Назвемо її

RK4m:

function (tout, yout| = RK4m(Zpfun, Mpfun, h, t, y)
% RK4m – Інтегрування ЗДР методом Рунге-Кута 4-го порядку.
% праві частини яких задані процедурою Zpfun і Mpfun.
% Вхідні змінні:
% Zpfun - рядок, що містить і'мя процедури обчислення правих частин ЗДР.
% звернення: z = fun(Mpfun, t, y), де Zpfun = 'fun',
% a Mpfun - рядок, що містить ім'я процедури, до якої
% звертається процедура fun;
% t - поточний момент часу
% у - вектор поточних значень змінних стану
% z - обчислені значення похідних z(i) = dy(i)/dt.
% h - крок інтегрування
% t - попередній момент часу
% у - попереднє значення вектора змінних стану.
% Вихідні змінні:
% tout - новий момент часу
% yout - обчислене значення вектора змінних стану через крок інтегрування
% Розрахунок проміжних значень похідних
kl = feval(Zpfun, Mpfun,t, у);
k2 = feval(Zpfun, Mpfun, t+h/3, y+h/3*kl);
k3 = feval(Zpfun, Mpfun, t+2*h/3, y+h*k2-h/3*kl);
k4 = feval(Zpfun, Mpfun, t+h, y+h*(k3+kl-k2));
% Розрахунок нових значень вектора змінних стану
tout = t + h;
yout = у + h*(kl + 3*k2 + 3*k3 + k4)/8;
% Кінець процедури RK4m
Ця форма подання процедури обчислення правих частин диференційних

рівнянь не є зручною, бо, по-перше, процедуру вигляду FM1 неможливо

використовувати у зверненні до процедур MatLAB ode23 та ode45 (останні