моделирование / Моделирование2
.docxМИНОБРНАУКИ РОССИИ
Федеральное государственное бюджетное образовательное учреждение высшего образования
«САРАТОВСКИЙ НАЦИОНАЛЬНЫЙ ИССЛЕДОВАТЕЛЬСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ
ИМЕНИ Н. Г. ЧЕРНЫШЕВСКОГО»
ОТЧЁТ ПО ЗАДАНИЮ №2
Студентка 4 курса 411 группы
Направления 02.03.02—Фундаментальная информатика и информационные технологии
Факультета КНиИТ
Филатова Ольга Владимировна
Проверил
ассистент И.А.Люкшин
Саратов2022
СОДЕРЖАНИЕ
Задание 3
Код программы и результаты выполнения 4
Задание
Задача 19. Построить фазовый портрет системы дифференциальных уравнений
Код программы и результаты выполнения
Программа написана на высокоуровневом интерпретируемом языке Octave.
Попробуем преобразовать эту систему уравнений 2 порядка в систему уравнений 1 порядка. Для этого обозначим x’ за a, а y' за b. Применим метод Эйлера с пересчетом для решения системы дифференциальных уравнений 1 порядка.
Код вспомогательной функции:
function [ t , resx , resy ] = eiler ( f1 , f2 , f3 , f4 , ts , N0)
t = ts (1) : 0.01 : ts (2) ;
h = 0.01;
a = zeros ([1 size(t ,2)]) ;
b = zeros ([1 size(t ,2)]) ;
a (1) = N0(1);
x(1) = N0(2);
b(1) = N0(1);
y(1) = N0(2);
for i = 1 : size(t ,2) - 1
ka1 = h * f1 (x(i ) , y(i ) , a( i ) , t ( i ) ) ;
kx1 = h * f2 (x(i ) , y(i ) , a( i ) , t ( i ) ) ;
kb1 = h * f3 (x(i ) , y(i ) , b(i ) , t ( i ) ) ;
ky1 = h * f4 (x(i ) , y(i ) , b(i ) , t ( i ) ) ;
ka2 = h * f1 (x(i ) + kx1 / 2, y(i ) , a( i ) + ka1 / 2, t ( i ) + h/2) ;
kx2 = h * f2 (x(i ) + kx1 / 2, y(i ) , a( i ) + ka1 / 2, t ( i ) + h/2) ;
kb2 = h * f3 (x(i ) , y(i ) + ky1 / 2, b(i ) + kb1 / 2, t ( i ) + h/2) ;
ky2 = h * f4 (x(i ) , y(i ) + ky1 / 2, b(i ) + kb1 / 2, t ( i ) + h/2) ;
ka3 = h * f1 (x(i ) + kx2 / 2, y(i ) , a( i ) + ka2 / 2, t ( i ) + h/2) ;
kx3 = h * f2 (x(i ) + kx2 / 2, y(i ) , a( i ) + ka2 / 2, t ( i ) + h/2) ;
kb3 = h * f3 (x(i ) , y(i ) + ky2 / 2, b(i ) + kb2 / 2, t ( i ) + h/2) ;
ky3 = h * f4 (x(i ) , y(i ) + ky2 / 2, b(i ) + kb2 / 2, t ( i ) + h/2) ;
ka4 = h * f1 (x(i ) + kx3, y(i ) , a( i ) + ka3, t ( i ) + h);
kx4 = h * f2 (x(i ) + kx3, y(i ) , a( i ) + ka3, t ( i ) + h);
kb4 = h * f3 (x(i ) , y(i ) + ky3, b(i ) + kb3, t ( i ) + h);
ky4 = h * f4 (x(i ) , y(i ) + ky3, b(i ) + kb3, t ( i ) + h);
a( i+1) = a( i ) + 1/6 * (ka1 + 2*ka2 + 2*ka3 + ka4);
x(i+1) = x(i ) + 1/6 * (kx1 + 2*kx2 + 2*kx3 + kx4);
b(i+1) = b(i ) + 1/6 * (kb1 + 2*kb2 + 2*kb3 + kb4);
y(i+1) = y(i ) + 1/6 * (ky1 + 2*ky2 + 2*ky3 + ky4);
endfor
resx = [x]' ;
resy = [y]' ;
endfunction
Код основной программы:
t0 = input("Введите t0:");
t1 = input("Введите t1:");
f1_2 = @(x,y,a, t )(-4 * x + sin (t));
f1_1 = @(x,y,a, t ) (a) ;
f2_2 = @(x,y,b,t )(4*x - 8*y);
f2_1 = @(x,y,b,t ) (b);
figure ;
hold on
[ t_gr , X_gr, Y_gr] = eiler (f1_1, f1_2, f2_1, f2_2, [t0 , t1 ], [0, 0]) ;
plot ( X_gr , Y_gr, 'r ' ) ;
xlabel ( ' x ' )
ylabel ( 'y ')
grid on
legend('x(t )', 'y(t )') ;
Скриншот результата выполнения программы: