close

Вход

Забыли?

вход по аккаунту

?

мой курсач 2

код для вставкиСкачать
Санкт-Петербургский государственный политехнический университет
факультет Технической кибернетики
кафедра "Системный анализ и управление"
Пояснительная записка
К курсовой работе
По дисциплине: "Вычислительная математика"
На тему: "Решение дифференциальных уравнений 3го порядка"
Выполнил:
Кладкевич М.В.
Проверил работу: А.Н. Кирсяев Санкт-Петербург
2012
Содержание:
1. Постановка задачи3
2. Исследование методов решения3
3. Аналитическое решение задачи4
4. Численное решение5
4.1. Исходные данные5
4.2. Описание функций программы5
4.3. Исходный код5
4.4. Отладка и тестирование программы7
4.5. Исследование14
5. Выводы14
1. Постановка задачи
Построить таблицу и график решения линейного дифференциального уравнения, заданнsv численным методом и аналитически, а также с использованием стандартной процедуры ode45. Дифференциальное уравнение:
a0 y'''(t) + a1 y''(t) + a2 y'(t) + y(t) = 10. Начальные условия: y(0) = y'(0) = y''(0) = 0. Интервал:
t Є [0, tm]. Сравнить полученные результаты с аналитическим решением. Исследовать зависимость погрешности решения от выбора шага интегрирования.
Исходные данные:
a0 = 2;
a1 = 1.02; a2 = 2.01; tm = 20 (с).
В ходе работы будет использоваться метод Адамса-Башфорта 3-го порядка.
2. Исследование методов решения
Общая формула Адамса-Башфорта: x_(k+1) = x_k+h∑_(r = 0)^m▒〖β_mr f_(k-r) 〗
Метод Адамса-Башфорта 3-го порядка является трехшаговым (поскольку использует информацию в трех точках), явным и имеет формулу: x_(k+1) = x_k+h/12∙(〖23f〗_k-〖16f〗_(k-1)+〖5f〗_(k-2) ).
Общая форма ЛДУ: a_0 y'''(t)+a_1 y''(t) +a_2 y'(t) +a_3 y(t) = u
ЛДУ уравнение можно привести к системе, а именно:
пусть y = x_1, y^' = x_1^' = x_2, y^'' = x_2^' = x_3, где x_1,x_2,x_3 - неизвестные числа;
подставим выше указанную замену в ЛДУ и выразим оттуда y^''' . Получим: y^''' = x_3^' = 1/a_0 (-a_3 x_1-a_2 x_2-a_1 x_3 )+u/a_0 ;
в результате получаем систему ЛДУ(СЛДУ): {█(x_1^' = x_2 @x_2^' = x_3@x_3^' = 1/a_0 (-a_3 x_1-a_2 x_2-a_1 x_3 )+u/a_0 )┤
Полученную СЛДУ можно переписать в матричном виде:
dx/dt = Ax+Bu, где A = ■(0&1&0@0&0&1@-a_3/a_0 &-a_2/a_0 &-a_1/a_0 ); B = ■(0@0@1/a_0 ); x = ■(x_1@x_2@x_3 )
Bu = f, следовательно можно записать для конкретного заданного уравнения:
f = ■(0@0@10/4);A = ■(0& 1&0@0& 0&1@-1/2&-2.01/2&-1.02/2); x_0 = ■(0@0@0);n = 3 m = dt/h-количество узлов.
Для расчёта дополнительных начальных условий требуется стартовый метод, например, одношаговый метод Рунге-Кутты.
xk+1 = xk+h*(k1+2*k2+2*k3+k4)/6, где
k1 = f(tk, xk)
k2 = f(tk+h/2, xk+h*k1/2)
k3 = f(tk+h/2, xk+h*k2/2)
k4 = f(tk+h, xk+h*k3)
h = (tk -t0 )/n
3. Аналитическое решение задачи 2y'''(t) + 1.02 y''(t) + 2.01 y'(t) + y(t) = 10
Из заданного дифференциального уравнения найдем корни характеристического многочлена.
С помощью следующей программы:
A = [0 1 0; 0 0 1; -1/2 -2.01/2 -1.02/2];
ei1 = eig(A)
Результат:
-0.5000 -0.0050 + 1.0000i
-0.0050 - 1.0000i
Так как собственные значения различаются больше, чем на порядок, то из этого следует, что эта система достаточно жесткая. Кроме того, ввиду того, что вещественные части всех собственных значений отрицательны, система асимптотически устойчива.
Соответствующие им частные решения ДУ:
Общее решение нашего ДУ:
Начальные условия: y(0) = y'(0) = y''(0) = 0. y'''(0) = 10 / a0 = 5
F.m
y = vpa(dsolve('2 * D3y + 1.02 * D2y + 2.01 * Dy + y = 10', 'y(0) = 0', 'D2y(0) = 0', 'Dy(0) = 0', 'x'))
Результат:
y =
10.0 - (1.9678714859437751004016064257028*cos(0.99998749992187402342224094390368*x))/exp(0.005*x) - (4.0259539391966028659938403448691*sin(0.99998749992187402342224094390368*x))/exp(0.005*x) - 8.0321285140562248995983935742972/exp(0.5*x)
4. Численное решение
4.1. Исходные данные
1) Дано:
a0 = 2;
a1 = 1.02; a2 = 2.01; tm = 20 (с);
y0 = [0;0;0].
2) Приведение ДУ к форме Коши 3) Расчет общего решения ДУ Y(k) аналитическим методом
4) Расчет с помощью ode45 5 )Расчет первых трёх значений с помощью метода Рунге-Кутты 6) Расчет с помощью метода Адамса-Башфорта 3 -го порядка
7) Вывод графиков расчета разными методами на экран
4.2. Описание функций программы
Основные функции и переменные, использованные в реализованной программе: function kursach2(a, y0, t0, tm, h) - функция, параметрами которой являются: коэффициенты в уравнении, начальные условия ЛДУ, интервал времени, шаг
X1 = R(t,a)- функция расчета методом Рунге-Кутты
[T,Y] = AB(tspan,a)- - функция, осуществляющая стартовый метод; выходными параметрами являются вектор-столбец Y и параметр T; входными параметрами являются начало и конец отрезка интегрирования; исходная матрица а
4.3. Исходный код
function kursach2(a, y0, t0, tm, h) a = [2, 1.02, 2.01]; % коэффициенты ДУ
A = [ 0, 1, 0; 0, 0, 1; - 1/a(1), - a(3)/a(1), - a(2)/a(1)]; % матрица A
t0 = 0;
tm = 20;
f = [0; 0; 10/a(1)]; h = 1; % шаг
y0 = [0; 0; 0]; % начальные условия n = (tm - t0)/h; % количество точек
for k = 1:n + 1 % цикл по всем точкам для задания аналитического решения
T(k) = t0 + (k - 1) * h; Y(k) = 10. -8.03212851405622489959839357430 * exp(-.5 * T(k)) - 4.02595393919660286599384034488 * exp(-.5e-2 * T(k)) * sin(.999987499921874023422240943905 * T(k)) - 1.96787148594377510040160642570 * exp( - .500000000000000000000000000000e-2 * T(k)) * cos(.999987499921874023422240943905 * T(k)); end
options = odeset('RelTol', 1e-10, 'AbsTol', 1e-5); % точность ode45
[T1, Y1] = ode45('F', T, y0, options); % решение ode45
% - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - [T2, Y2] = AB(T, a); % метод Адамса - Башфорта для шага h
[T3, Y3] = AB(T/2, a); % метод Адамса - Башфорта для шага h/2
% построение графиков
figure(1); subplot(2, 1, 1); plot(T, (Y1(:, 1)));
grid on;
title('Решение ode45');
subplot(2, 1, 2); plot(T, (Y' - Y2(:, 1)));
grid on;
title('Погрешность ode45');
figure(2); subplot(2, 1, 1); plot(T, (Y' - Y1(:, 1)));
grid on;
title('Метод Адамсa-Башвортa');
subplot(2, 1, 2); plot(T, (Y2(:, 1)));
grid on;
title('Погрешность Адамса-Башворта');
figure(3);
plot(T, Y); grid on;
title('Аналитическое');
function X1 = R(t, a) % Рунге
A = [ 0, 1, 0; 0, 0, 1; - 1/a(1), - a(3)/a(1), - a(2)/a(1)];
f = [0; 0; 10/a(1)];
y = [0; 0; 0]; % начальные условия
y1 = A * y + f;
X1(:, 1) = y; h = t; % for k = 1:2 % цикл
k1 = A * y + f; % k2 = A * (y + h/2 * k1) + f; k3 = A * (y + h/2 * k2) + f; k4 = A * (y + h * k3) + f; y = y + h * (k1 + 2 * k2 + 2 * k3 + k4)/6; y1 = A * y + f; X1(:, k + 1) = y; end % конец % явный метод адамса - башфорта
function [T, Y] = AB(tspan, a) % диапазон времени
n = length(tspan); % кол - во шагов и точек разбиения
h = tspan(n)/(n - 1); % шаг
% X(:, 1) = y0'; % начальные значения
f_own = (R(h, a)); % рассчитываем с помощью метода Рунге - Кутты первые три значения
X(:, 1) = f_own(:, 1); X(:, 2) = f_own(:, 2); X(:, 3) = f_own(:, 3); % записываем их в столбцы
for k = 3:n - 1
X(:, k + 1) = X(:, k) + h/12 * (23 * F(tspan(k), X(:, k)) - 16 * F(tspan(k - 1), X(:, k - 1)) + 5 * F(tspan(k - 2), X(:, k - 2)) ); % метод адамса - башфорта end
T = tspan; Y = X'; 4.4. Отладка и тестирование программы
Для сравнения погрешностей при разных шагах и анализа устойчивости, будут проводиться вычисления с различным шагом. Найдём h критическое для метода Адамса-Башфорта 3го порядка:
h*| λ | < 6/11
λ = -0.5
h < 6/(11*| λ |)
h критическое ~1,1
h критическое получилось больше единицы, поэтому критическим шагом в данной работе будем считать h = 1.
Ниже представлены графики для различных h и анализ погрешностей.
h = 0,005
h = 0,01
h = 0,05 h = 0,1
h = 0,5 h = 1
4.5. Исследование
Шаг hПогрешность метода А-БПогрешность ode450,0053*10-62*10-50,012,5*10-52*10-50,053*10-32*10-50,10,032*10-50,52,52*10-511,3*1042*10-5
После проведения тестирования видно, что погрешность метода Адамса-Башфорта увеличивается с увеличением ширины шага. На шаге h < 0,01 она становится меньше погрешности стандартного ode45.
Можно сделать вывод, что при равномерном выборе величины шага метод Адамса-Башфорта является недостаточно устойчивым для большой величины шага. Точности решения можно добиться, если подбирать величину шага с течением времени, обращая внимание на требования устойчивости.
5. Выводы
В ходе работы была реализована программа решения ДУ методом Адамса-Башфорта 3-го порядка, использующая также метод Рунге-Кутты.
Была проведена отладка и тестирование программы, где полученное с помощью данного метода решение было сравнено с решением с помощью встроенного ode45. В результате работы была найдена зависимость погрешности интегрирования от выбранного шага.
Документ
Категория
Рефераты
Просмотров
91
Размер файла
457 Кб
Теги
мой, курсач
1/--страниц
Пожаловаться на содержимое документа