[Link].
ru
4.2. Численные методы решение краевой задачи для ОДУ
Примером краевой задачи является двухточечная краевая задача для
обыкновенного дифференциального уравнения второго порядка.
y ' ' = f ( x, y , y ' ) (4.28)
с граничными условиями, заданными на концах отрезка [a, b] .
y (a) = y 0
(4.29)
y (b) = y1
Следует найти такое решение y (x) на этом отрезке, которое принимает на
концах отрезка значения y 0 , y1 . Если функция f ( x, y, y ' ) линейна по аргументам
y , y ' , то задача (4.28),(4.29) - линейная краевая задача, в противном случае –
нелинейная.
Кроме граничных условий (4.29) называемых граничными условиями
первого рода, используются еще условия на производные от решения на концах -
граничные условия второго рода:
y ' (a ) = yˆ 0
(4.30)
y ' (b) = yˆ1
или линейная комбинация решений и производных - граничные условия третьего
рода:
αy (a) + β y ' (a) = yˆ 0
, (4.31)
δy (b) + γy ' (b) = yˆ1
где α , β , δ , γ - такие числа, что α + β ≠ 0, δ + γ ≠ 0 .
Возможно на разных концах отрезка использовать условия различных
типов.
В данном пособии рассматриваются два приближенных метода решения краевой
задачи:
- метод стрельбы (пристрелки);
- конечно-разностный метод.
4.2.1. Метод стрельбы
1
[Link]
Суть метода заключена в многократном решении задачи Коши для
приближенного нахождения решения краевой задачи.
Пусть надо решить краевую задачу (4.28), (4.29) на отрезке [a, b] . Вместо
исходной задачи формулируется задача Коши с уравнением (4.28) и с
начальными условиями
y (a) = y 0
, (4.32)
y ' (b) = η
где η - некоторое значение тангенса угла наклона касательной к решению в точке
x = a.
Положим сначала некоторое начальное значение параметру η = η 0 , после чего
решим каким либо методом задачу Коши (4.28),(4.32). Пусть y = y 0 ( x, y 0 ,η 0 )
решение этой задачи на интервале [a, b] , тогда сравнивая значение функции
y 0 (b, y 0 ,η 0 ) со значением y1 в правом конце отрезка можно получить
информацию для корректировки угла наклона касательной к решению в левом
конце отрезка. Решая задачу Коши для нового значения η = η1 , получим другое
решение со значением y1 (b, y 0 ,η1 ) на правом конце. Таким образом, значение
решения на правом конце y (b, y 0 ,η ) будет являться функцией одной переменной
η. Задачу можно сформулировать таким образом: требуется найти такое
значение переменной η * , чтобы решение y (b, y 0 ,η * ) в правом конце отрезка
совпало со значением y1 из (4.29). Другими словами решение исходной задачи
эквивалентно нахождению корня уравнения
Φ (η ) = 0 , (4.33)
где Φ (η ) = y (b, y 0 ,η ) − y1 .
Уравнение (4.33) является “алгоритмическим” уравнением, так как левая
часть его задается с помощью алгоритма численного решения соответствующей
задачи Коши. Но методы решения уравнения (4.33) аналогичны методам
решения нелинейных уравнений, изложенным в разделе 2. Следует заметить, что
так как невозможно вычислить производную функции Φ (η ) , то вместо метода
2
[Link]
Ньютона следует использовать метод секущих, в котором производная от функции
заменена ее разностным аналогом. Данный разностный аналог легко вычисляется
по двум приближениям, например η k и η k +1 . Следующее значение искомого
корня определяется по соотношению
η j +1 − η j
η j + 2 = η j +1 − Φ (η j +1 ) (4.34)
Φ (η j +1 ) − Φ (η j )
Итерации по формуле (4.34) выполняются до удовлетворения заданной
точности.
Пример 4.9. Методом стрельбы решить краевую задачу y ' ' = e x + sin y с
граничными условиями 1-го рода y (0) = 1, y (1) = 2 на отрезке [0,1] .
Решение
Заменой переменных z = y ' сведем дифференциальное уравнение второго
порядка к системе двух дифференциальных уравнений первого порядка.
⎧ y' = z
⎨
⎩ z ' = e + sin y
x
Задачу Коши для системы с начальными условиями на левом конце
y (0) = 1, y ' (0) = η будем решать методом Рунге-Кутта 4-го порядка точности с
шагом h = 0.1 до удовлетворения условия на правом конце
y (1.0,1.0,η k ) − 2.0 = Φ(η k ) ≤ ε , где ε = 0.0001, и y (1.0,1.0,η k ) - значение решения
задачи Коши в правом конце отрезка при b = 1.0, y (0) = y 0 = 1.0, η k - значение
первой производной к решению в левом конце отрезка на k – ой итерации.
Примем в качестве первых двух значений параметра η следующие: η 0 =1.0,
η1 =0.8. Дважды решим задачу Коши с этими параметрами методом Рунге-Кутта с
шагом h =0.1, получим два решения y (1.0,1.0,η 0 ) = 3.168894836, y(1.0,1.0,η1 ) =
2.97483325. Вычислим новое приближение параметра η по формуле (4.34)
3
[Link]
0 . 8 − 1 .0
η 2 = 0 .8 − (2.97483325 − 2.0) = −0.204663797 ;
2.97483325 − 3.168894836
Решая задачу Коши с параметром η 2 , получим решение y (1.0,1.0,η 2 ) =
1.953759449 и так далее.
− 0.204663797 − 0.8
η 3 = −0.204663797 − (1.953759449 − 2.0) = −0.159166393 ;
1.953759449 − 2.97483325
y (1.0,1.0,η 3 ) = 2.001790565; Φ(η 3 ) = 0.001790565 ≥ ε ;
− 0.159166393 − ( −0.204663797 )
η 4 = −0.159166393 − ( 2.001790565 − 2.0) = −0.160862503 ;
2.001790565 − 1.953759449
y (1.0,1.0,η 4 ) = 2.000003115; Φ(η 4 ) = 0.000003115 ≤ ε ;
Вычисления заносим в таблицу 4.15
Таблица 4. 15
j ηj y (1.0,1.0,η j ) Φ (η j )
0 +1.000000000 3.168894836 1.168894836
1 +0.800000000 2.974483325 0.974483325
2 -0.204663797 1.953759449 0.046240551
3 -0.159166393 2.001790565 0.001790565
4 -0.160862503 2.000003115 0.000003115
Приближенным решением краевой задачи будем считать табличную
функцию, полученную в результате решения задачи Коши с параметром η 4 и
приведенную в таблице 4.16.
Таблица 4.16
x k 0. 0.100 0.200 0.300 0.400 0.500 0.600 0.700 0.800 0.900 1.00
0 00 00 00 00 00 00 00 00 00 0
y k 1.0 0.993 1.006 1.039 1.094 1.1743 1.279 1.4123 1.5752 1.770 2.00
28 01 42 97 4 44 6 8 45 0
4
[Link]
4.2.2. Конечно-разностный метод решения краевой задачи
Рассмотрим двухточечную краевую задачу для линейного
дифференциального уравнения второго порядка на отрезке [a, b]
y ' '+ p ( x ) y '+ q ( x ) y = f ( x ) (4.35)
y (a ) = y 0 , y (b) = y1 . (4.36)
Введем разностную сетку на отрезке [a, b] Ω ( h ) = {x k = x 0 + hk } , k = 0,1,...., N ,
h = b − a / N . Решение задачи (4.35),(4.36) будем искать в виде сеточной функции
y ( h ) = {y k , k = 0,1,...., N }, предлагая, что решение существует и единственно. Введем
разностную аппроксимацию производных следующим образом:
y k +1 − y k −1
y k' = + O(h 2 ) ;
2h
y k +1 − 2 y k + y k −1
y k'' = + O(h 2 ) ; (4.37)
h2
Подставляя аппроксимации производных из (4.37) в (4.35),(4.36) получим
систему уравнений для нахождения y k :
⎧ y0 = ya
⎪ y − 2y + y y − y k −1
⎪ k +1
⎨ 2
k k −1
+ p( x k ) k +1 + q( x k ) y k = f ( x k ), k = 1, N − 1 (4.38)
⎪ h 2h
⎪⎩ y N = y b
Приводя подобные и учитывая, что при задании граничных условий
первого рода два неизвестных y 0 , y N уже фактически определены, получим
систему линейных алгебраических уравнений с трехдиагональной матрицей
коэффициентов
5
[Link]
⎧ p ( x1 )h p ( x1 )h
⎪( −2 + h 2
q ( x 1 ) y 1 + (1 + ) y 2 = h 2
f ( x1 ) − ( 1 − ) ya
2 2
⎪
⎪ p( x k )h p ( x k )h
⎨(1 − ) y k −1 + (−2 + h 2 q ( x k )) y k + (1 + ) y k +1 = h 2 f ( x k ) ,k=2,...,N-2
⎪ 2 2
⎪ p ( x N −1 )h p ( x N −1 )h
⎪(1 − ) y N −1 + (−2 + h 2 q ( x N −1 )) y N −1 = h 2 f ( x N −1 ) − (1 + ) yb
⎩ 2 2
(4.39)
Для системы (4.39) при достаточно малых шагах сетки h и q ( x k ) < 0
выполнены условия преобладания диагональных элементов
p( x k )h p( x k )h
− 2 + h 2 q( x k ) > 1 − + 1+ , (4.39)
2 2
что гарантирует устойчивость счета и корректность применения метода прогонки
для решения этой системы.
В случае использования граничных условий второго и третьего рода
аппроксимация производных проводится с помощью односторонних разностей
первого и второго порядков.
y1 − y 0
y 0' = + O ( h) ;
h
y N − y N −1
y N' = + O ( h) (4.40)
h
− 3 y 0 + 4 y1 − y 2
y 0' = + O(h 2 ) ;
2h
y N − 2 − 4 y N −1 + 3 y N
y N' = + O(h 2 ) ; (4.41)
2h
В случае использования формул (4.40) линейная алгебраическая система
аппроксимирует дифференциальную задачу в целом только с первым порядком
(из-за аппроксимации в граничных точках), однако сохраняется трех
диагональная структура матрицы коэффициентов. В случае использования
формул (4.41) второй порядок аппроксимации сохраняется везде, но матрица
линейной системы не трехдиагональная.
6
[Link]
⎧ y ′′ − xy ′ − y = 0
⎪
Пример 4.10. Решить краевую задачу ⎨ y (0) = 1 с шагом h = 0.2.
⎪ y ′(1) + 2 y (1) = 0
⎩
Здесь p ( x) = x , q( x) = 1 , f ( x) = 0 , N = 5, x0 = 0, x1 = 0.2, x 2 = 0.4, x3 = 0.6,
x 4 = 0.8, x5 = 1.0
Во всех внутренних узлах отрезка [0,1] после замены производных их разностными
аналогами получим
(1 − 0.1x k ) y k −1 + (−2.04) y k + (1 + 0.1x k ) y k +1 = 0 , k = 1,...., 4
На левой границе y 0 = 1 , на правой границе аппроксимируем производную
односторонней разностью 1-го порядка:
y5 − y 4
+ 2 y5 = 0 .
0.2
С помощью группировки слагаемых, приведения подобных членов и
подстановки значений x k и с учетом y 0 = 1 получим систему линейных
алгебраических уравнений.
⎧− 2.04 y1 + 1.02 y 2 = −0.98
⎪0.96 y − 2.04 y + 1.04 y = 0
⎪⎪ 1 2 3
⎨0.94 y 2 − 2.04 y 3 + 1.06 y 4 = 0
⎪0.92 y − 2.04 y + 1.08 y = 0
⎪ 3 4 5
⎪⎩+ y 4 − 1.4 y 5 = 0
В данной трехдиагональной системе выполнено условие преобладания
диагональных элементов и можно использовать метод прогонки (раздел 1.1.2).
В результате решения системы методом прогонки получим следующие
значения: y5 = 0.2233205 , y 4 = 0.31265 , y 3 = 0.43111 , y 2 = 0.58303 , y1 = 0.77191.
Решением краевой задачи является табличная функция
Таблица 4.17
k 0 1 2 3 4 5
xk 0 0.2 0.4 0.6 0.8 1.0
yk 1.0 0.77191 0.58303 0.43111 0.31265 0.22332
7
[Link]
Найдите больше информации на сайте Учитесь.ру ([Link])!