Исследование неявного метода Эйлера для линейной системы ОДУ с постоянным и переменным шагом (стр. 2 из 2)
t2e=t2(1:max(size(t2))-1);
t3e=t3(1:max(size(t3))-1);
plot(t1e,e1,t2e,e2,t3e,e3)
pause
disp('Решаем постоянный шаг:')
pause;
% Постоянный шаг
[tc1,yc1,ec1]=nrk1('a1','b','u',0,3.5,[0.2;0.2],0.1);
[tc2,yc2,ec2]=nrk1('a1','b','u',0,3.5,[0.2;0.2],0.01);
[tc3,yc3,ec3]=nrk1('a1','b','u',0,3.5,[0.2;0.2],0.005);
plot(tc1,yc1,tc2,yc2,tc3,yc3)
pause
t1ec=tc1(1:max(size(tc1))-1);
t2ec=tc2(1:max(size(tc2))-1);
t3ec=tc3(1:max(size(tc3))-1);
plot(t1ec,ec1,t2ec,ec2,t3ec,ec3)
pause
[tc1,yc1,ec1]=nrk1('a1','b','u',0,3.5,[0.1;0.1],0.1);
[tc2,yc2,ec2]=nrk1('a1','b','u',0,3.5,[0.5;0.5],0.1);
[tc3,yc3,ec3]=nrk1('a1','b','u',0,3.5,[1;1],0.1);
plot(tc1,yc1,tc2,yc2,tc3,yc3)
pause
t1ec=tc1(1:max(size(tc1))-1);
t2ec=tc2(1:max(size(tc2))-1);
t3ec=tc3(1:max(size(tc3))-1);
plot(t1ec,ec1,t2ec,ec2,t3ec,ec3)
pause
- Нежесткая матрица А:
function A=a();
A=[-5/6 1/3;
1/3 -1/3];
Жесткая матрица А:
function A=a();
A=[-50 50;
50 -50.1];
Нежесткая матрица А:
function B=b();
B=[5/2; 0];
Матрица U:
function U=u(t);
U=[1];
Переменные
- Выходные переменные:
tout – выходной вектор времени;
yout – выходное значение функций (соответствует вектору времени);
eout – выходное значение ошибок определения функций (высчитано по методу Рунге).
- Внутренние переменные:
A, B, U – соответствующие матрицы (вектора) из формулы Х = АХ+ВU(t).
t – текущее время;
y – текущие значения функций;
ym1 – предыдущие значения функций;
yp1 – следующие значения функций;
yt –значения функций, высчитанное с половинчатым шагом (для определения ошибки по методу Рунге);
h – текущий шаг;
h1 – предыдущий шаг;
n – вектор с единицами с размерностями как y;
I – единичная матрица;
eam – текущая ошибка, вычисляемая для определения допустимого шага;
- Входы:
funA, funB, funU – внешние функции вычисления A, B и U;
t0, tfinal – начальное и конечное значение времени;
y0 – начальное значение функций (с которых начинаем решение системы ОДУ);
ep – допустимая ошибка;
trace – выводить или нет на экран промежуточные значения;
h – шаг в методе с постоянным шагом.
Описание алгоритма:Сначала производятся все начальные инициализации и установки. Затем запускается основной цикл расчета значений функций. В нем для переменного шага высчитывается критерий. Затем по этому критерию определяется сам шаг. Цикл повторяется. Для постоянного шага производится непосредственно пересчет значения функции. Также, для определения ошибки по формуле Рунге, производится вычисление функции через два половинных шага. На выходе получается значения времени, функций и ошибок. "Длина" массива ошибок на 1 меньше, чем "длины" других массивов (в силу специфики его вычисления).
5. Описание тестовых задач
В программе тестирования сначала находится решение нежесткой системы с разными начальными значениями y, шагами и способами: с переменным и постоянным шагами. Строятся совокупные графики функций и ошибок для каждого случая. Затем решается нежесткая система с разными допустимыми ошибками и разными начальными значениями y и разными шагами. Строятся совокупные графики функций и ошибок для каждого случая. Для способа с переменным шагом изменяем только начальные значения, для способа с постоянным шагом – начальные значения, и величины шага.
График1:
График2.
График 3.
График 4.
График5.
График 6.
График 7.
График 8.
График 9.
График 10.
эйлер линейный программа интегрирование
График 11.
График 12.
6. Анализ результатов. Выводы
Ниже показаны полученные графики функций и ошибок:
График1:
Var, XX0=[0.1;0.1] X0=[0.5;0.5] X0=[1;1]
График2:
Var, EX0=[0.1;0.1] X0=[0.5;0.5] X0=[1;1]
График3:
Const, X h=0.1, X0=[0.2;0.2] h=0.01, X0=[0.2;0.2] h=0.005, X0=[0.2;0.2]
График4:
Const, E h=0.1, X0=[0.2;0.2] h=0.01, X0=[0.2;0.2] h=0.005, X0=[0.2;0.2]
График5:
Const, X h=0.1, X0=[0.1;0.1] h=0.1, X0=[0.5;0.5] h=0.1, X0=[1;1]
График6:
Const, E h=0.1, X0=[0.1;0.1] h=0.1, X0=[0.5;0.5] h=0.1, X0=[1;1]
Пояснения:
Var– метод с переменным шагом;
Const – метод с постоянным шагом.
X0=[0.1;0.1], X0=[0.5;0.5], X0=[1;1] отражают цветовую гамму на графике.
X– графики функций;
E – графики ошибок.
h – шаг, X0 – начальное значение X.
Рисунки с 7 -12 аналогичны рисункам 1-6 (соответственно) по начальным данным и отличаются лишь жесткостью решаемой системы.
Выводы
Видим, что при увеличении начального значения происходит лишь параллельный сдвиг решения, однако график ошибки перемещается в противоположную сторону (при увеличении начального значения она уменьшается и наоборот), причем непропорционально. Следовательно, при увеличении начального значения X ошибка уменьшается, однако с увеличением t это изменение становится все менее и менее заметным. Это характерно как для метода с переменным шагом, так и для метода с постоянным шагом. Кроме того, для метода с постоянным шагом увеличение шага приводит к уменьшению точности решения. Что касается жесткой системы, то для нее характерно прямолинейное конечное решение. С уменьшением шага ошибка, посчитанная по методу Рунге, уменьшается. При увеличении начального X происходит параллельный сдвиг графика конечного решения. Если в случае нежесткой системы графики решений разных составляющих X находятся далеко друг от друга, то при решении жесткой системы получаются очень близкие графики.
Заключение
В ходе выполнения работы был изучен неявный метод Эйлера для решения линейных систем ОДУ. В ходе реализации проекта и проведения тестирования была проверена справедливость теоретических выкладок. Получены сведения о зависимости точности интегрирования от величины и способа выбора шага, а также от начальных значений переменных.
Используемая литература
1. Сарычева О.М. Численные методы. – Новосибирск, 1995г. – 65с.
2. Бахвалов Н.С. Численные методы. Ч1.- М: Наука, 1975г. – 632с., илл.
3. Копченова Н.В., Марон И.А. Вычислительная математика в примерах и задачах . – М: Наука, 1972г. - 368с