<<

Краткие теоретические сведения

Точное (аналитическое) решение (общее или частное) дифференциального уравнения подразумевает получение искомого решения (функции y(x)) в виде выражения от элементарных функций.
Это возможно далеко не всегда даже для уравнений первого порядка.

Численное решение ДУ (частное) заключается в вычислении функции y(x) и ее производных в некоторых заданных точках , лежащих на определенном отрезке. То есть, фактически, решение ДУ n-го порядка вида получается в виде следующей таблицы чисел (столбец значений старшей производной вычисляется подстановкой значений в уравнение):

X y y' y(n-1)
x1 y(x1) y'(x1) y(n-1)(x1)
x2 y(x2) y'(x2) y(n-1)(x2)
xN y(xN) y'(xN) y(n-1)(xN)

Например, для дифференциального уравнения первого порядка таблица решения будет представлять собой два столбца – x и y.

Множество значений абсцисс в которых определяется значение функции, называют сеткой, на которой определена функция y(x). Сами координаты при этом называют узлами сетки. Чаще всего, для удобства, используются равномерные сетки, в которых разница между соседними узлами постоянна и называется шагом сетки или шагом интегрирования дифференциального уравнения

или , i = 1, …, N

Для определения частного решения необходимо задать дополнительные условия, которые позволят вычислить константы интегрирования.

Причем таких условий должно быть ровно n. Для уравнений первого порядка – одно, для второго - 2 и т.д. В зависимости от способа их задания при решении дифференциальных уравнений существуют три типа задач:

· Задача Коши (начальная задача): Необходимо найти такое частное решение дифференциального уравнения, которое удовлетворяет определенным начальными условиям, заданным в одной точке:

то есть, задано определенное значение независимой переменной (х0), и значение функции и всех ее производных вплоть до порядка (n-1) в этой точке. Эта точка (х0) называется начальной. Например, если решается ДУ 1-го порядка, то начальные условия выражаются в виде пары чисел (x0, y0)

Такого рода задача встречается при решении ОДУ, которые описывают, например, кинетику химических реакций. В этом случае известны концентрации веществ в начальный момент времени (t = 0), и необходимо найти концентрации веществ через некоторый промежуток времени (t). В качестве примера можно так же привести задачу о теплопереносе или массопереносе (диффузии),

Численные методы решения задачи Коши ОДУ первого порядка

Рассмотрим некоторые численные методы решения задачи Коши (начальной задачи) обыкновенных дифференциальных уравнений первого порядка. Запишем данное уравнение в общем виде, разрешенном относительно производной (правая часть уравнения не зависит от первой производной):

(1)

Необходимо найти значения функции y в заданных точках сетки , если известны начальные значения , где есть значение функции y(x) в начальной точке x0.

Преобразуем уравнение умножением на dx

И проинтегрируем левую и правую части между i-ым и i+1-ым узлами сетки.

(2)

Мы получили выражение для построения решения в i+1 узле интегрирования через значения x и y в i-ом узле сетки. Сложность, однако, заключается в том, что интеграл в правой части есть интеграл от неявно заданной функции, нахождение которого в аналитическом виде в общем случае невозможно. Численные методы решения ОДУ различным способом аппроксимируют (приближают) значение этого интеграла для построения формул численного интегрирования ОДУ.

Из множества разработанных для решения ОДУ первого порядка методов рассмотрим метод Эйлера, Он достаточно прост и даёт начальное представление о подходах к решению данной задачи в рамках численного решения задачи Коши.

Метод Эйлера

Исторически первым и наиболее простым способом численного решения задачи Коши для ОДУ первого порядка является метод Эйлера. В его основе лежит аппроксимация производной отношением конечных приращений зависимой (y) и независимой (x) переменных между узлами равномерной сетки:

где yi+1 это искомое значение функции в точке xi+1.

Если теперь преобразовать это уравнение, и учесть равномерность сетки интегрирования, то получится итерационная формула, по которой можно вычислить yi+1 , если известно yi в точке хi:

(3)

Сравнивая формулу Эйлера с общим выражением, полученным ранее, видно, что для приближенного вычисления интеграла в (2) в методе Эйлера используется простейшая формула интегрирования - формула прямоугольников по левому краю отрезка.

Графическая интерпретация метода Эйлера также не представляет затруднений (см. рисунок ниже). Действительно, исходя из вида решаемого уравнения (1) следует, что значение есть значение производной функции y(x) в точке x=xi - , и, таким образом, равно тангенсу угла наклона каcательной, проведенной к графику функции y(x) в точке x=xi.

Из прямоугольного треугольника на рисунке можно найти

,

откуда и получается формула Эйлера. Таким образом, суть метода Эйлера заключается в замене функции y(x) на отрезке интегрирования прямой линией, касательной к графику в точке x=xi. Если искомая функция сильно отличается от линейной на отрезке интегрирования, то погрешность вычисления будет значительной. Ошибка метода Эйлера прямо пропорциональна шагу интегрирования:

Ошибка ~ h

Процесс вычислений строится следующим образом. При заданных начальных условиях x0 и y0 можно вычислить

Таким образом, строится таблица значений функции y(x) с определенным шагом (h) по x на отрезке [x0, xN]. Ошибка в определении значения y(xi) при этом будет тем меньше, чем меньше выбрана длина шага h (что определяется точностью формулы интегрирования).

При больших h метод Эйлера весьма неточен. Он дает все более точное приближение при уменьшении шага интегрирования. Если отрезок [xi, xi+1] слишком велик, то каждый участок [xi, xi+1] разбивается на N отрезков интегрирования и к каждому их них применяется формула Эйлера с шагом , то есть шаг интегрирования h берется меньше шага сетки, на которой определяется решение.

Пример1:

Используя метод Эйлера, построить приближенное решение для следующей задачи Коши:

На сетке с шагом 0,1 в интервале [0, 1] (4)

Решение:

Данное уравнение уже записано в стандартном виде, резрешенном относительно производной искомой функции.

Поэтому, для решаемого уравнения имеем

Примем шаг интегрирования равным шагу сетки h = 0,1.

При этом для каждого узла сетки будет вычислено только одно значение (N=1). Для первых четырех узлов сетки вычисления будут следующими:

И т. д.

Полные результаты (с точностью до пятого знака после запятой) приведены в таблице 1 в третьей колонке - h=0,1 (N=1). Во второй колонке таблицы для сравнения приведены значения, вычисленные по аналитическому решению данного уравнения .

Во второй части таблицы приведена относительная погрешность полученных решений. Видно, что при h=0,1 погрешность весьма велика, достигая 100% для первого узла x=0,1.

Таблица 1 Решение уравнения (4) методом Эйлера (для колонок указан шаг интегрирования и число отрезков интегрирования N между узлами сетки)

bgcolor="#99CCFF">0,0015625
x Точное решение 0,1 0,05 0,025 0,00625 0,0007813 0,0001953
1 2 4 16 64 128 512
0 0,000000 0,000000 0,000000 0,000000 0,000000 0,000000 0,000000 0,000000
0,1 0,004837 0,000000 0,002500 0,003688 0,004554 0,004767 0,004802 0,004829
0,2 0,018731 0,010000 0,014506 0,016652 0,018217 0,018603 0,018667 0,018715
0,3 0,040818 0,029000 0,035092 0,037998 0,040121 0,040644 0,040731 0,040797
0,4 0,070320 0,056100 0,063420 0,066920 0,069479 0,070110 0,070215 0,070294
0,5 0,106531 0,090490 0,098737 0,102688 0,105580 0,106294 0,106412 0,106501
0,6 0,148812 0,131441 0,140360 0,144642 0,147779 0,148554 0,148683 0,148779
0,7 0,196585 0,178297 0,187675 0,192186 0,195496 0,196314 0,196449 0,196551
0,8 0,249329 0,230467 0,240127 0,244783 0,248202 0,249048 0,249188 0,249294
0,9 0,306570 0,287420 0,297214 0,301945 0,305423 0,306284 0,306427 0,306534
1 0,367879 0,348678 0,358486 0,363232 0,366727 0,367592 0,367736 0,367844
Относительные погрешности вычисленных значений функции при различных h
x h 0,1 0,05 0,025 0,00625 0,0015625 0,0007813 0,0001953
N 1 2 4 16 64 128 512
0,1 100,00% 48,32% 23,76% 5,87% 1,46% 0,73% 0,18%
0,2 46,61% 22,55% 11,10% 2,74% 0,68% 0,34% 0,09%
0,3 28,95% 14,03% 6,91% 1,71% 0,43% 0,21% 0,05%
0,4 20,22% 9,81% 4,83% 1,20% 0,30% 0,15% 0,04%
0,5 15,06% 7,32% 3,61% 0,89% 0,22% 0,11% 0,03%
0,6 11,67% 5,68% 2,80% 0,69% 0,17% 0,09% 0,02%
0,7 9,30% 4,53% 2,24% 0,55% 0,14% 0,07% 0,02%
0,8 7,57% 3,69% 1,82% 0,45% 0,11% 0,06% 0,01%
0,9 6,25% 3,05% 1,51% 0,37% 0,09% 0,05% 0,01%
1 5,22% 2,55% 1,26% 0,31% 0,08% 0,04% 0,01%

Уменьшим шаг интегрирования вдвое, h = 0.05, в этом случае для каждого узла сетки вычисление будет проводиться за два шага (N=2). Так, для первого узла x=0,1 получим:

И так далее, до конца отрезка.

Из таблицы 1 (четвертая колонка, N=2) видно, что погрешность решения резко снизилась, примерно вдвое, хотя и осталась по-прежнему, значительной.

При шаге интегрирования h=0,025 для каждого узла сетки необходимо выполнить 4 вычисления по формуле Эйлера в промежуточных точках (N=4).

(Для других узлов значения приведены в таблице 1, колонка N=4)

В таблице 1 приведены для сравнения вычисления для некоторых других значений N, вплоть до 512. Видно, что точность решения возрастает весьма медленно при уменьшении шага интегрирования, необходимо брать очень маленький шаг для достижения приемлемой точности (и, следовательно, много раз вычислять значение F(x,y)). Поэтому метод Эйлера практически не используется в вычислительной практике.

Усовершенствованный метод Эйлера ( метод Эйлера с пересчетом).

Точность метода Эйлера можно повысить, если воспользоваться для аппроксимации интеграла более точной формулой интегрирования – формулой трапеций.

(5)

Данная формула оказывается неявной относительно yi+1 (это значение есть и в левой и в правой части выражения), то есть является уравением относительно yi+1, решать которое можно, например, численно, применяя какой-либо итерационный метод (в таком виде его можно рассматривать как итерационную формула метода простой итерации). Однако, можно поступить иначи и приблизительно вычислить значение функции в узле i+1 с помощью обычной формулы Эйлера:

,

которое затем использовать при вычислении по (5).

Таким образом получается метод Эйлера с пересчетом. Для каждого узла интегрирования производится следующая цепочка вычислений

(6)

Благодаря более точной формуле интегрирования, погрешность метода Эйлера с пересчетом пропорциональна уже квадрату шага интегрирования.

Ошибка ~ h2

Пример2

Проведем вычисления для уравнения (4) с помощью метода Эйлера с пересчетом .

При шаге интегрирования h=0,1 в первом узле сетки x1 получим:

Что намного точнее значения, полученного методом Эйлера при том же шаге интегрирования. В таблице 2 ниже приведены сравнительные результаты вычислений при h= 0,1 методов Эйлера и Эйлера с пересчетом

Таблица 2 Решение уравнения методами Эйлера и Эйлера с пересчетом

x Точное метод Эйлера с пересчетом Метод Эйлера
y отн. погрешность y отн. погрешность
0 0,000000 0,00000 0,00000
0,1 0,004837 0,00500 3,36% 0,00000 100,00%
0,2 0,018731 0,01903 1,57% 0,01000 46,61%
0,3 0,040818 0,04122 0,98% 0,02900 28,95%
0,4 0,070320 0,07080 0,69% 0,05610 20,22%
0,5 0,106531 0,10708 0,51% 0,09049 15,06%
0,6 0,148812 0,14940 0,40% 0,13144 11,67%
0,7 0,196585 0,19721 0,32% 0,17830 9,30%
0,8 0,249329 0,24998 0,26% 0,23047 7,57%
0,9 0,306570 0,30723 0,21% 0,28742 6,25%
1 0,367879 0,36854 0,18% 0,34868 5,22%

Отметим существенное увеличение точности вычислений метод Эйлера с пересчетом по сравнению с методом Эйлера. Так, для узла x=0,1 относительное отклонение значения функции, определенного метод Эйлера с пересчетом , оказывается в 30 (!) раз меньше. Такая же точность вычислений по формуле Эйлера достигается при числе отрезков интегрирования N примерно 30. Следовательно, при использовании метода Эйлера с пересчетом при одинаковой точности вычислений понадобится примерно в 15 раз меньше времени ЭВМ, чем при использовании метода Эйлера.

Проверка устойчивости решения

Решение ОДУ в некоторой точке xi называется устойчивым, если найденное в этой точке значение функции yi мало изменяется при уменьшении шага интегрирования. Для проверки устойчивости, таким образом, надо провести два расчета значения (yi) – с шагом интегрирования h и при уменьшенной (например, двое) величине шага

В качестве критерия устойчивости можно использовать малость относительного изменения полученного решения при уменьшении шага интегрирования (ε – наперед заданная малая величина)

Такая проверка может осуществляться и для всех решений на всем интервале значений x. Если условие не выполняется, то шаг снова делится пополам и находится новое решение и т.д. до получения устойчивого решения.

Контрольные вопросы

1. Метод Эйлера

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

Литература

1. Данко П.Е., Попов А.Г., Кожевникова Т.Я. «Высшая математика в упражнениях и задачах», части 1, 2. М. ОНИКС 21 век, Мир и Образование. 2002г.
2. Омельченко В.П., Курбатова Э.В. «Математика». Ростов-на-Дону. Феникс. 2007г.
3. Яковлев Г.Н. «Алгебра и начала анализа», часть 2. М. Наука. 88г.

Литература

1.
Башмаков М.И. «Алгебра и начала анализа 10-11класс». М. Дрофа. 2000г.

2. Валуце И.И., Дилигул Г.Д. «Математика для техникумов». М. Наука. 80г.

3. Гмурман В.Е. "Теория вероятностей и математическая статистика". М. Высшая школа. 99г.

4. Данко П.Е., Попов А.Г., Кожевникова Т.Я. «Высшая математика в упражнениях и задачах», части 1, 2. М. ОНИКС 21 век, Мир и Образование. 2002г.

5. Омельченко В.П., Курбатова Э.В. «Математика». Ростов-на-Дону. Феникс. 2007г.

6. Филимонова Е.В. «Математика». Ростов-на-Дону. Феникс. 2005г.

7. Шипачёв В.С. «Основы высшей математики». М. Высшая школа. 94г.

8. Шипачёв В.С. «Сборник задач по высшей математике». М. Высшая школа. 94г.

9. Яковлев Г.Н. «Алгебра и начала анализа», часть 1. М. Наука. 87г.

10. Яковлев Г.Н. «Алгебра и начала анализа», часть 2. М. Наука. 88г.

<< |
Источник: Лабгаева Эмма Владимировна. Методические указания для студентов по проведению практических занятий по дисциплине «Математика». 2007

Еще по теме Краткие теоретические сведения:

  1. ИЗБРАННЫЕ ТЕОРЕТИЧЕСКИЕ ВОПРОСЫ СУДЕБНОЙ МЕДИЦИНЫ
  2. Теоретические основы русской лексикографии
  3. §1. Разработка теоретических основ и особенности развития правового регулирования общественных отношений в условиях НЭПа
  4. Раздел  II. ПРАВО (Общая теория права. Право: общетеоретические понятияи определения)
  5. 6.2. Краткие теоретические сведения.
  6. С О Д Е Р Ж А Н И Е
  7. Введение
  8. Краткие теоретические сведения
  9. Краткие теоретические сведения
  10. Краткие теоретические сведения
  11. Краткие теоретические сведения
  12. Краткиетеоретические сведения
  13. Краткие теоретические сведения
  14. Краткие теоретические сведения
  15. Краткие теоретические сведения
  16. Краткие теоретические сведения