Численные методы решения систем дифференциальных уравнений

1

СОДЕРЖАНИЕ

1. Решение систем линейных уравнений

а) Схема Халецкого

б) Метод Зейделя и условия сходимости

2. Методы решения нелинейных уравнений

а) Метод хорд

б) Метод Ньютона (метод касательных)

в) Метод итерации

3. Приближенное вычисление обыкновенных дифференциальных уравнений

а) Метод Рунге-Кутта

б) Экстраполяционные методы Адамса

в) Метод Милна

г) Краевые задачи для обыкновенных дифференциальных уравнений

4. Приближенные методы решения дифференциальных уравнений с частными производными

а) Классификация дифференциальных уравнений второго порядка

б) Постановка краевых задач

в) Метод конечных разностей (метод сеток)

г) Разностные схемы для решения уравнения теплопроводности

д) Разностные схемы для решения уравнения колебания струны

Список литературы

1. Решение систем линейных уравнений

Системы линейных уравнений (СЛУ) имеют в вычислениях очень большое значение, так как к ним может быть приведено приближенное решение широкого круга задач. Так, основными источниками возникновения СЛУ являются теория электрических цепей, уравнения балансов и сохранения в механике, гидравлике и т.д. Существует несколько способов решения таких систем, которые в основном делятся на два типа: 1) точные методы, представляющие собой конечные алгоритмы для вычисления корней системы, 2) итерационные методы, позволяющие получать корни системы с заданной точностью путем сходящихся бесконечных процессов. Заметим, что даже результаты точных методов являются приближенными из-за неизбежных округлений. Для итерационных процессов также добавляется погрешность метода.

Пример системы линейных уравнений:

Или в матричном виде: ,

где

матрица коэффициентов системы;

- вектор неизвестных; - вектор свободных членов.

Схема Халецкого

Запишем систему линейных уравнений в матричном виде: ,

где A=[a_ij] - квадратная матрица порядка n и

, - векторы-столбцы.

Представим матрицу Aв виде произведения нижней треугольной матрицы B=[b_ij] и верхней треугольной матрицы C=[c_ij] с единичной диагональю , где

и .

Тогда элементы b_ij и c_ij определяются по формулам

и

Отсюда искомый вектор x может быть вычислен из уравнений и .

Так как матрицы B и C - треугольные, то системы легко решаются:

и

Из этих двух формул видно, что числа y_i выгодно вычислять вместе с коэффициентами c_ij. Этот метод получил название схемы Халецкого. В схеме применяется обычный контроль с помощью сумм. Если матрица A - симметрическая

a_ij=a_ji, то

Пример. Решить систему

Решение.

В первый раздел таблицы впишем матрицу коэффициентов системы, ее свободные члены и контрольные суммы. Далее так как , то первый столбец из раздела 1 переносится в первый столбец раздела II. Чтобы получить первую строку раздела II, делим все элементы первой строки раздела I на элемент, в нашем случае на 3.

Имеем:

; ; ; ; .

Переходим к заполнению второго столбца раздела II, начиная со второй строки. Пользуясь формулами, определяем

:

;

;

.

Далее определяя по формулам, заполняем вторую сетку для раздела II:

Затем переходим к третьему столбцу, вычисляя его элементы и по формулам и т.д., пока не будет заполнена вся таблица раздела II. Таким образом, заполнение раздела II происходит способом “елочки”: столбец - строка, столбец - строка и т.д.

В разделе Ш, пользуясь формулами, определяем и .

Текущий контроль осуществляется с помощью столбца ?, над которым производятся те же действия, что и над столбцом свободных членов.

I							3	1	-1	2	6	11
I							-5	1	3	-4	-12	-17
I							2	0	1	1	1	3
I							1	-5	3	3	3	-1
II							3	0.333333	-0.333333	2	2	3.666667
II							-5	2.666667	-0.25	0.25	-0.75	0.5
II							2	-0.666667	2	-1.25	-1.75	-2
II							1	-5.333333	6	2.5	3	4
III											2	1
III											-0.75	-1
III											-1.75	2
III											3	3

Метод Зейделя и условия сходимости

Этот метод представляет собой модификацию метода простой итерации. Его смысл заключается в том, что при вычислении (k+1)-го приближения неизвестной x_i учитываются уже вычисленные ранее (k+1)-е приближения x_1, x₂, ..., x_i-1. Пусть дана приведенная линейная система (i = 1, 2, …n). Выберем произвольно начальные приближения корней , стараясь, чтобы они в какой-то мере соответствовали искомым неизвестным x₁, x₂, x₃, ..., x_n.

Предположим, что k-е приближение корней известно, тогда в соответствии с идеей метода будем строить (k+1)-е приближение по следующим формулам:

Обычно процесс Зейделя сходится быстрее, чем метод простой итерации. Бывает, что процесс Зейделя сходится, когда простая итерация расходится и т.п. Правда, бывает и наоборот. Во всяком случае, достаточные условия сходимости для метода простой итерации достаточны и для сходимости метода Зейделя. То есть процесс итерации сходится, если выполнено одно из условий

1) или 2) .

Пример. Методом Зейделя решить систему уравнений

Решение. Приведем эту систему к виду, удобному для итерации,

В качестве нулевых приближений корней возьмем:

; ; .

Применяя процесс Зейделя, последовательно получим:

и т.д.

Результаты вычислений с точностью до четырех знаков помещены в таблице:

0	1,2000	0,0000	0,0000
1	1,2000	1 ,0600	0,9480
2	0,9992	1,0054	0,9991
3	0,9996	1.0001	1,0001
4	1 ,0000	1,0000	1,0000
5	1 ,0000	1,0000	1,0000

Точные значения корней: .

2. Методы решения нелинейных уравнений

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

Метод хорд

Пусть дано уравнение f(x)=0, где функция f(x) определена и непрерывна на интервале

[a, b] и f(a)f(b)<0. Пусть для определенности f(a)<0 и f(b)>0. Разделим отрезок [a, b] в отношении - f(a):f(b). Это даст нам приближенное значение корня x₁ = a + h₁, где

.

Далее этот прием применяем к одному из отрезков [a, x₁] или [x₁, b], на концах которого функция f(x) имеет противоположные знаки. Аналогично находим второе приближение x₂ и т.д. Геометрически этот способ эквивалентен замене кривой y = f(x) хордой, проходящей через точки А[a, f(a)] и B[b, f(b)].



f(b)

о

f(a)

0 x

о x₃ x₂ x₁ b=x₀ a=x₀ x₁ x₂ b

0 f(a) x

a

f(b)

Действительно, уравнение хорды АВ имеет вид

При х = х₁ и y = 0, получим

Полагая, что на отрезке [a, b] вторая производная f''(x) сохраняет постоянный знак, метод хорд сводится к двум различным вариантам.

Из рис. 1 видно, что конец а неподвижен и последовательные приближения: x₀=b;

образуют ограниченную монотонно убывающую последовательность, причем a<о<…<x_n+1<x_n<…<x₁<x₀.

Из рис. 2 видно, что неподвижен конец b и последовательные приближения: x₀=a;

образуют ограниченную монотонно возрастающую последовательность, причем

x₀<x₁<x₂<…<x_n<x_n+1<…<о<b.

Таким образом, для вычисления корня уравнения имеем две различные вычислительные формулы. За неподвижный конец выбираем тот конец, для которого знак функции f(x) совпадает со знаком второй производной f''(x).

Пример. Найти положительный корень уравнения

с точностью до 0,002.

Решение. Прежде всего отделяем корень. Так как и , то искомый корень лежит в интервале . Полученный интервал велик, поэтому разделим его пополам. Так как то .Последовательно применяя формулы, будем иметь:

Так как

и при

имеем

то можно принять:

Таким образом, , где . Заметим, что точный корень уравнения есть .

Метод Ньютона (метод касательных)

Пусть корень о уравнения f(x)=0, отделен на отрезке [a, b], причем первая и вторая производные f'(x) и f''(x) непрерывны и сохраняют определенные знаки при . Найдя какое-нибудь n-ое приближение корня , мы можем уточнить его по методу Ньютона следующим образом. Пусть о=x_n+h_n, где h_n - величина малая. Отсюда по формуле Тейлора получим: f(x_n + h_n) ? f(x_n)+h_n fў(x_n)=0.

Следовательно,

.

Подставив полученное выражение в формулу о=x_n+h_n, найдем следующее значение корня:



Графическое нахождение корня методом Ньютона (рис. 3).

Если в качестве начального приближения выбрать точку х₀=В₀ , то процесс быстро сходится. Если же выбрать точку х₀=А, то х₁[a, b], и процесс нахождения корня расходится. Рекомендуется: в качестве х₀ выбирать такую точку, где f(x₀)f''(x₀)>0.

Пример. Вычислить методом Ньютона отрицательный корень уравнения

с пятью верными знаками.

Решение. Полагая в левой части уравнение получим . Следовательно, искомый корень находится в интервале . Сузим найденный интервал. Так как то . В этом последнем интервале и .Так как и , то можем принять за начальное приближение . Последовательные приближения вычисляем по следующей схеме:

Останавливаясь на , проверяем знак значения . Так как , то, и любое из этих чисел дает искомое приближение.

0	-11	3453	-5183	0,7
1	-10,3	134,3	-4234	0,03
2	-10,27	37,8	-4196	0,009
3	-10,261	0,2	-	-

Метод итерации

Заменим уравнение f(x)=0 эквивалентным x=ц(x). Выберем некоторое начальное приближение x₀ и вычислим дальнейшие приближения по формулам x₁= ц(x₀), x₂= ц(x₁), …, x_n= ц(x_n-1). Если последовательность x_n имеет предел, то итерационный процесс

x_n= ц(x_n-1) (n=1, 2, …) называется сходящимся. Пусть функция ц(x) непрерывна. Переходя к пределу в равенстве x_n= ц(x_n-1), получим

Следовательно, является корнем уравнения x=ц(x) и может быть вычислен по формуле x_n= ц(x_n-1) (n=1, 2, …) с любой точностью. Для данного метода существуют две теоремы:

Теорема 1. Пусть корень уравнения x=ц(x), а также его последовательные приближения x₀, x_n= ц(x_n-1) (n=1, 2, …) содержатся в интервале [a, b] и на [a, b]. Тогда справедливы утверждения:

1. итерационный процесс x_n= ц(x_n-1) сходится к корню уравнения ;

_2. ;

_3. .

Следствие 1. Если требуется найти корень с точностью е, то кончаем итерационный процесс тогда, когда

<е, т.е. когда .

Следствие 2. Так как =ц() и =ц(x_n-1), то -x_n= ц()-ц(x_n-1). По теореме Лагранжа . Из этого следует, что если ц'(x)>0 на (a, b), то последовательные приближения x_n= ц(x_n-1) (n=1, 2, …) сходятся к корню монотонно; если ц'(x)<0 на (a, b), то последовательные приближения колеблются около корня.

Теорема 2. Если на [a, b], а корень и начальное приближение x₀ находятся на более узком отрезке [б, в], где , то справедливы заключения теоремы 1.

Привести уравнение f(x)=0 к виду x=ц(x) таким образом, чтобы получить сходящийся итерационный процесс, можно различными способами. Рассмотрим два из них:

1) уравнение f(x)=0 равносильно при л?0 уравнению лf(x)=0 и уравнению x= лf(x)+x. Обозначим лf(x)+x через ц(x), получим x= ц(x). Параметр л подберем так, чтобы функция ц'(x)= лf'(x)+1 на [a, b] была по модулю меньше единицы.

2) если , то итерационный процесс расходится. Заменим уравнение x=ц(x) эквивалентным ему уравнением x=ш(x), где ш(x) - функция, обратная функции ц(x). Так как , то итерационный процесс x_n=ш(x_n-1) будет сходящимся.

Пример. Методом итерации найти корень уравнения 5x-8lnx=8 с точностью 0,01.

Решение. Запишем уравнение в виде и построим соответствующие графики:

Уравнение имеет два корня: . За начальные приближения возьмем z₀=0,5 и x₀=3,5. Для уточнения запишем . Здесь

Следовательно, итерационный процесс сходится. Погрешность оценим по формуле , результаты вычислений приведены в таблице:

n	x	1+lnx
0 1 2 3 4	3,5 3,605 3,651 3,672 3,682	2,253 2,282 2,295 2,301	3,605 3,651 3,672 3,682	------ 0,105 0,046 0,021 0,010

Так как ц'(z₀)?3>1, то итерационный процесс расходится. Найдем функцию , обратную функции ц(x). Так как , то итерационный процесс будет сходится. , результаты вычислений приведены в таблице:

n	zn
0 1 2	0,5 0,503 0,503	-0,688 -0,686 ------	0,503 0,504 ------	------ 0,0015 0,0005

3.Приближенное вычисление обыкновенных дифференциальных уравнений

Обыкновенным дифференциальным уравнением называется равенство , в котором - независимая переменная, изменяющаяся в некотором отрезке , а - неизвестная функция от , которую и надо найти. Различают два типа обыкновенных дифференциальных уравнений - уравнения без начальных условий и уравнения с начальными условиями. Уравнения без начальных условий - это как раз то, что было только что определено. А уравнение с начальными условиями - это записанное выше уравнение относительно функции , но в котором требуется найти лишь такую функцию , которая удовлетворяет при некотором следующим условиям:

,

т.е. в точке функция и ее первые производных принимают наперед заданные значения. В этой ситуации число называется порядком уравнения.

5. Приближенное вычисление обыкновенных дифференциальных уравнений

Обыкновенным дифференциальным уравнением называется равенство , в котором - независимая переменная, изменяющаяся в некотором отрезке , а - неизвестная функция от , которую и надо найти. Различают два типа обыкновенных дифференциальных уравнений - уравнения без начальных условий и уравнения с начальными условиями. Уравнения без начальных условий - это как раз то, что было только что определено. А уравнение с начальными условиями - это записанное выше уравнение относительно функции , но в котором требуется найти лишь такую функцию , которая удовлетворяет при некотором следующим условиям:

,

т.е. в точке функция и ее первые производных принимают наперед заданные значения. В этой ситуации число называется порядком уравнения.

Метод Рунге-Кутта

Изложим идею метода на примере:

Интегрируя это уравнение в пределах от x до x + h (0 < h <1), получим равенство

которое посредством последнего интеграла связывает значения решения рассматриваемого уравнения в двух точках, удаленных друг от друга на расстояние шага h. Для удобства записи данного выражения используем обозначение ?y=y(x+h)-y(x) и замену переменной интегрирования t=x+ah. Окончательно получим:



Указав эффективный метод приближенного вычисления интеграла в выражении

,

мы получим при этом одно из правил численного интегрирования уравнения

Постараемся составить линейную комбинацию величин j_i, i = 0, 1, ..., q, которая будет являться аналогом квадратурной суммы и позволит вычислить приближенное значение приращения Dy:

где

Метод четвертого порядка для q = 3, имеет вид

где

Особо широко известно другое вычислительное правило Рунге-Кутта четвертого порядка точности:

где

Метод Рунге-Кутта имеет погрешность четвертого порядка (~ h⁴ ).

Правило Рунге. Если приближенный метод имеет порядок погрешности m, то погрешность можно приближенно оценить по формуле

В формуле O(x_i) - главный член погрешности, и - приближенные решения в точке x_i, найденные с шагом h и 2h соответственно.

Экстраполяционные методы Адамса

Широко распространенным семейством многошаговых методов являются методы Адамса. Простейший из них, получающийся при , совпадает с рассмотренным ранее методом Эйлера первого порядка точности. В практических расчетах чаще всего используется вариант метода Адамса, имеющий четвертый порядок точности и использующий на каждом шаге результаты предыдущих четырех. Именно его и называют обычно методом Адамса. Рассмотрим этот метод.

Пусть найдены значения в четырех последовательных узлах . При этом имеются также вычисленные ранее значения правой части . В качестве интерполяционного многочлена можно взять многочлен Ньютона. В случае постоянного шага конечные разности

для правой части в узле имеют вид

.

Тогда разностная схема четвертого порядка метода Адамса запишется в виде

.

Сравнивая метод Адамса с методом Рунге -- Кутта той же точности, отмечаем его экономичность, поскольку он требует вычисления лишь одного значения правой части на каждом шаге. Но метод Адамса неудобен тем, что невозможно начать счет по одному лишь известному значению . Расчет может быть начат лишь с узла . Значения необходимые для вычисления , нужно получить каким-либо другим способом , что существенно усложняет алгоритм. Кроме того, метод Адамса не позволяет изменить шаг в процессе счета; этого недостатка лишены одношаговые методы.

Метод Милна

Пусть на отрезке [a, b] требуется найти численное решение дифференциального уравнения с начальным условием . Разобьем отрезок [a, b] на n равных частей точками , где h=(b-a)/n - шаг интегрирования. Используя начальные данные, находим каким-либо способом последовательные значения искомой функции y(x). Таким образом, становится известным . Приближения и для следующих значений последовательно находятся по формулам Милна

- где .

Абсолютная погрешность значения приближенно равна

.

Пример. Дано дифференциальное уравнение y'=y-x, удовлетворяющие начальному условию x₀=0, y(x₀)=1,5. Вычислить с точность до 0,01 значение решения этого уравнения при x=1,5.

Решение. Выберем начальный шаг вычисления. Из условия h⁴<0,01 получим h=0,25 Составим таблицу

i	xi	yi	y'i=f(xi, yi)=yi-xi		y'i= f(xi, yi)=yi-xi		еi
0 1 2 3 4 5 6	0 0,25 0,50 0,75 1,00 1,25 1,50	1,5000 1,8920 2,3243 2,8084	1,5000 1,6420 1,8243 2,0584	3,3588 3,9947 4,7402	2,3588 2,7447 3,2402	3,3590 3,9950 4,7406	7*10-5 10-5 1,4*10-5

Получаем ответ y=(1,5)=4,74.

Краевые задачи для обыкновенных дифференциальных уравнений

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

, а дополнительные условия ставятся более, чем в одной точке отрезка интегрирования уравнений (в этом случае порядок системы не может быть меньше второго): , , .

Общая классификация методов решения краевых задач: существуют точные, приближенные и численные методы.

4. Приближенные методы решения дифференциальных уравнений с частными производными

Кроме обычных дифференциальных уравнений существуют так называемые дифференциальные уравнения с частными производными. Далее они будут рассмотрены более подробно.

Классификация дифференциальных уравнений второго порядка

Рассмотрим уравнение второго порядка

,

где - функции и . Говорят, что указанное уравнение в области принадлежит гиперболическому типу, если в этой области . Если , то уравнение в области принадлежит параболическому типу. Если , то уравнение принадлежит эллиптическому типу.

Уравнение

называется каноническим уравнением гиперболического типа.

Уравнение

называется каноническим уравнением параболического типа.

Уравнение

называется каноническим уравнением эллиптического типа.

Дифференциальное уравнение

называется уравнением характеристик уравнения

.

Если последнее уравнение гиперболического типа, то уравнение характеристик имеет два интеграла: т.е. существуют два семейства вещественных характеристик.

С помощью замены переменных , дифференциальное уравнение

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

.

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

В этом случае осуществляем замену переменных , , где -- какая-нибудь функция, для которой

.

После замены переменных получаем уравнение

.

Для уравнения эллиптического типа интегралы уравнения характеристик имеют вид

, где и -- вещественные функции.

Полагая и , приводим уравнение

к виду

.

Постановка краевых задач

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

1. Задача Коши для бесконечной области. Рассмотрим эту задачу на примере

уравнения колебания струны и уравнения теплопроводности.

Рассмотрим процесс колебания тонкой бесконечной струны под действием непрерывно распределенной внешней силы с плотностью f. Предположим, что сила действует в одной плоскости - плоскости колебания струны (x, u), а струна является гибкой упругой нитью. Пусть величина натяжения, возникающая в струне вследствие ее изгиба, подчиняется закону Гука, а сами колебания достаточно малы. Тогда величина смещения u (x, t) удовлетворяет уравнению колебания струны: . Для однозначности процесса необходимо задать еще начальное смещение и начальное распределение скоростей. Математически это соответствует заданию начальных условий: . Требуется найти классическое решение уравнения, удовлетворяющие начальным условиям. Сформулированная таким образом задача называется задачей Коши для гиперболического уравнения.

Исследуем теперь процесс распределения температуры в тонком бесконечном стержне. Предполагается, что тепловой поток подчиняется закону Фурье, а изменение температуры тела пропорционально количеству теплоты, сообщаемой телу. Предположим, что внутри стержня может выделяться и поглощаться теплота, характеризуемая плотностью тепловых источников f. Тогда распределение температуры в стержне описывается уравнением теплопроводности: . Для однозначного задания процесса необходимо указать начальное распределение температуры. Это соответствует заданию начального условия: . Требуется найти классическое решение уравнения, удовлетворяющие начальным условиям. Сформулированная таким образом задача называется задачей Коши для параболического уравнения.

2. Стационарная задача (задача без начальных данных). Рассмотрим установившийся

режим распределения температуры в ограниченной тонкой пластине произвольной формы с гладкой границей. Пусть функция u(x, y) выражает температуру каждой точки пластины. При обычных законах распространения тепла функция u(x, y) удовлетворяет уравнению Пуассона: , где функция а задает плотность тепловых источников пластины. В случае отсутствия источника (f=0) данное уравнение называется уравнением Лапласа: . Для однозначного описания процесса необходимо задать тепловой режим на границе пластины. Это может быть сделано с помощью задания распределения температуры на границе или распределения теплового потока. Возможен также режим теплового равновесия излучающего тела с окружающей средой.

В зависимости от теплового режима на границе получаются три граничных условия для функции u(x, y). Пусть Г - граница рассматриваемой области D - определения уравнения Лапласа. Математическая формулировка граничных условий может быть задана в следующем виде:

граничное условие I рода: ;

граничное условие II рода: ;

граничное условие III рода: .

Производная берется по внешней нормали к кривой Г; л>0 - коэффициент теплопроводности; ц₀, ц₁, ц₂ - заданные на Г функции, причем ц₂ есть произведение коэффициента теплопроводности на температуру внешней среды, соприкасающейся с телом.

Таким образом, краевая задача заключается в том, чтобы найти классическое решение уравнения Пуассона или Лапласа, удовлетворяющее одному из граничных условий.

3. Смешанная краевая задача. Рассмотрим задачу распространения тепла в тонком

стержне единичной длины. Поместим один из концов в точку x=0, а другой - в точку x=1. Распределение температуры в таком стержне в течение некоторого интервала времени 0<t<T описывается уравнением , с начальным условием , а для единственности решения в этом случае необходимо еще задать температурный режим на концах стержня. Это можно сделать с помощью граничных условий, аналогичных тем, которые были сформулированы для уравнений Пуассона и Лапласа.

Граничное условие I рода (на конце стержня x=0 заданна температура): .

Граничное условие II рода (на конце стержня x=0 задан тепловой поток): .

Граничное условие III рода: .

Для другого конца стержня x=1 правые части граничных условий заменяются соответственно на ш₀(t), ш₁(t), ш₂(t). Заметим, что начальное и граничное условия должны удовлетворять так называемым условиям сопряжения, т.е. при условии I рода u₀(0)=ц₀(0), при условии II рода u_0x(0)=ц₁(0), при условии III рода -u_0x(0)+лu₀(0)=ц₂(0). Аналогичные условия сопряжения должны выполнятся и на другом конце стержня x=1.

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

Метод конечных разностей (метод сеток)

Численные методы, основанные на разностной аппроксимации производных называется разностным методом, методом конечных разностей или методом сеток.

Пусть заданно линейное дифференциальное уравнение, записанное в символическом виде: . Здесь u - искомое решение уравнения; L - некоторый дифференциальный оператор, сокращенно обозначающий соответствующую дифференциальную операцию; f - правая часть уравнения (заданная функция).

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

Разностный метод решения этих двух задач можно представить в виде двух этапов:

1. построение разностной схемы, аппроксимирующей данную непрерывную задачу;

2. получение решения разностной задачи и оценка погрешности этого решения.

Для построения разностной схемы первым шагом является замена области непрерывного изменения аргументов областью дискретного их изменения - сеточной областью , т.е.множеством точек (x_n, y_m), называемых узлами сетки. Для квадрата сеточную область можно построить следующим образом. Проведем прямые . Множество точек пересечения этих прямых и составит сеточную область, а сами точки образуют узлы сетки. Всякая функция , определенная на ссеке , называется сеточной функцией и обозначается .

Второй шаг в построении разностной схемы состоит в аппроксимации дифференциального выражения Lu некоторым разностным выражением, а функцию непрерывного аргумента f - сеточной функцией, т.е. в построение некоторого разностного аналога для данного уравнения, при данных краевых условиях.

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

Где L_h и ц_h - разностные операторы, аппроксимирующие соответственно L и l; х_h - искомая сеточная функция, аппроксимирующая решение u; f_h, ц_h - заданные сеточные функции, аппроксимирующие f и ц.

Совокупность разносных уравнений, аппроксимирующих исходную задачу - есть разностная схема. Рассмотрим их подробнее на примерах уравнения теплопроводности и колебания струны.

Разностные схемы для решения уравнения теплопроводности (параболический тип)

Рассмотрим первую краевую задачу для уравнения теплопроводности в прямоугольнике . Требуется найти непрерывное в решение задачи:

В области введем прямоугольную равномерную сетку {x_n, t_k} с шагом h=1/N по координате x и с шагом ф=T/M по координате t:

.

Производные левой части уравнения

аппроксимируем следующим разностными выражениями:

В соответствии с данной аппроксимацией построим два разностных аналога уравнения с неизвестной сеточной функцией х_hф:

Здесь - значение некоторой сеточной функции f_hф, соответствующей правой части уравнения . Для первой разностной схемы

,

а для второй -

.

Начальное и граничное условия для первой краевой задачи аппроксимируются точно:

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

Полагая r=ф/h² получим

- для первой разностной схемы,

- - для второй разностной схемы.

Анализ показывает, что погрешность аппроксимации схем есть .

Разностные схемы для решения уравнения колебания струны (гиперболический тип)

Рассмотрим первую краевую задачу для уравнения колебания струны в прямоугольнике . Требуется найти непрерывное в решение задачи:

Применение метода конечных разностей к решению задачи по существу мало чем отличается от его применения к уравнению теплопроводности. Область покрывается сеткой . Отличие заключается в приближении второй производной по переменной t:

.

Разностная аппроксимация принимает вид

.

Начальные условия аппроксимируются следующим образом:

.

Граничные условия аппроксимируются точно так же, как и для уравнения теплопроводности:

.

Значение является фиктивным неизвестным, которое можно определить по формуле:

, где г=ф/h.

Анализ показывает, что погрешность аппроксимации схем есть

.

Список литературы

1. Демидович Б.П., Марон И.А. Основы вычислительной математики. Наука, 1970.

2. Минкова Р.М., Вайсбурд Р.А. Методы вычислительной математики. УПИ, 1981.

3. Боглаев Ю.П. Вычислительная математика и программирование. Высшая школа, 1990.

4. Кацман Ю.Я. Прикладная математика. Численные методы. ТПУ, 2000.