Методы решения краевых задач, в том числе «жестких» краевых задач
Методы Алексея Юрьевича Виноградова
1 Введение
На примере системы дифференциальных уравнений цилиндрической оболочки ракеты – системы обыкновенных дифференциальных уравнений 8-го порядка (после разделения частных производных).
Система линейных обыкновенных дифференциальных уравнений имеет вид:
Y 
(x) = A(x) ∙ Y(x) + F(x),
где Y(x) – искомая вектор-функция задачи размерности 8х1, Y (x) – производная искомой вектор-функции размерности 8х1, A(x) – квадратная матрица коэффициентов дифференциального уравнения размерности 8х8, F(x) – вектор-функция внешнего воздействия на систему размерности 8х1.
Здесь и далее вектора обозначаем жирным шрифтом вместо черточек над буквами
Краевые условия имеют вид:
U∙Y(0) = u,
V∙Y(1) = v,
где
Y(0) – значение искомой вектор-функции на левом крае х=0 размерности 8х1, U – прямоугольная горизонтальная матрица коэффициентов краевых условий левого края размерности 4х8, u – вектор внешних воздействий на левый край размерности 4х1,
Y(1) – значение искомой вектор-функции на правом крае х=1 размерности 8х1, V – прямоугольная горизонтальная матрица коэффициентов краевых условий правого края размерности 4х8, v – вектор внешних воздействий на правый край размерности 4х1.
В случае, когда система дифференциальных уравнений имеет матрицу с постоянными коэффициентами A=const, решение задачи Коши имеет вид [Гантмахер]:
Y(x) = e 
∙ Y(x 
) + e 
∙ $IMAGE6$ e $IMAGE7$∙ F(t) dt,
где
e = E + A(x-x ) + A $IMAGE10$ (x-x ) $IMAGE10$/2! + A $IMAGE13$ (x-x ) $IMAGE13$/3! + …,
где E это единичная матрица.
Матричная экспонента ещё может называться матрицей Коши или матрициантом и может обозначаться в виде:
K(x←x ) = K(x - x ) = e .
Тогда решение задачи Коши может быть записано в виде:
Y(x) = K(x←x ) ∙ Y(x ) + Y*(x←x ) ,
где Y*(x←x ) = e ∙ $IMAGE6$ e $IMAGE7$∙ F(t) dt это вектор частного решения неоднородной системы дифференциальных уравнений.
2 Случай переменных коэффициентов
Этот вариант рассмотрения переменных коэффициентов проверялся в кандидатской диссертации.
Из теории матриц [Гантмахер] известно свойство перемножаемости матричных экспонент (матриц Коши):
e $IMAGE26$= e $IMAGE27$∙ e $IMAGE28$ ∙ … ∙ e $IMAGE29$ ∙ e $IMAGE30$,
K(x $IMAGE31$←x ) = K(x $IMAGE31$←x $IMAGE34$) ∙ K(x $IMAGE34$←x $IMAGE36$) ∙ … ∙ K(x $IMAGE37$←x $IMAGE38$) ∙ K(x $IMAGE38$←x ).
В случае, когда система дифференциальных уравнений имеет матрицу с переменными коэффициентами A=A(x), решение задачи Коши предлагается искать при помощи свойства перемножаемости матриц Коши. То есть интервал интегрирования разбивается на малые участки и на малых участках матрицы Коши приближенно вычисляются по формуле для постоянной матрицы в экспоненте. А затем матрицы Коши, вычисленные на малых участках, перемножаются:
K(x $IMAGE31$←x ) = K(x $IMAGE43$←x $IMAGE34$) ∙ K(x $IMAGE34$←x $IMAGE36$) ∙ … ∙ K(x $IMAGE37$←x $IMAGE38$) ∙ K(x $IMAGE38$←x ),
где матрицы Коши приближенно вычисляются по формуле:
K(x $IMAGE51$←x $IMAGE52$) = e $IMAGE53$, где ∆x $IMAGE43$= x $IMAGE55$- x $IMAGE31$.
3 Формула для вычисления вектора частного решения неоднородной системы дифференциальных уравнений
Эта очень простая формула еще не обсчитана на компьютерах. Вместо неё обсчитывалась значительно ранее выведенная и гораздо более сложная формула, приведенная в:
Численный метод переноса краевых условий для жестких дифференциальных уравнений строительной механики Журнал "ММ", Том: 14 (2002), Номер: 9, 3 стр. 1409-003r.pdf
Вместо формулы для вычисления вектора частного решения неоднородной системы дифференциальных уравнений в виде [Гантмахер]:
Y*(x←x ) = e ∙ $IMAGE6$ e $IMAGE7$∙ F(t) dt
предлагается использовать следующую формулу для каждого отдельного участка интервала интегрирования и тогда вектор частного решения на всем интервале будет складываться из векторов, вычисленных по формуле:
Y*(x $IMAGE61$←x $IMAGE31$) = Y*(x $IMAGE61$- x $IMAGE31$) = K(x $IMAGE61$- x $IMAGE31$) ∙ $IMAGE67$K(x $IMAGE31$- t) ∙ F(t) dt .
Правильность приведенной формулы подтверждается следующим:
Y*(x $IMAGE61$- x $IMAGE31$) = e $IMAGE71$∙ $IMAGE67$e $IMAGE73$∙ F(t) dt ,
Y*(x $IMAGE61$- x $IMAGE31$) = $IMAGE67$e $IMAGE71$∙e $IMAGE73$∙ F(t) dt ,
Y*(x $IMAGE61$- x $IMAGE31$) = $IMAGE67$e $IMAGE82$∙ F(t) dt ,
Y*(x $IMAGE61$- x $IMAGE31$) = $IMAGE67$e $IMAGE86$∙ F(t) dt ,
Y*(x $IMAGE61$- x $IMAGE31$) = e $IMAGE89$∙ $IMAGE67$e $IMAGE91$∙ F(t) dt ,
Y*(x←x $IMAGE31$) = e ∙ $IMAGE94$ e $IMAGE7$∙ F(t) dt,
что и требовалось подтвердить.
Вычисление вектора частного решения системы дифференциальных уравнений производиться при помощи представления матрицы Коши под знаком интеграла в виде ряда и интегрирования этого ряда поэлементно:
Y*(x $IMAGE61$←x $IMAGE31$) = Y*(x $IMAGE61$- x $IMAGE31$) = K(x $IMAGE61$- x $IMAGE31$) ∙ $IMAGE67$K(x $IMAGE31$- t) ∙ F(t) dt =
= K(x $IMAGE61$- x $IMAGE31$) ∙ $IMAGE67$ (E + A(x $IMAGE31$- t) + A $IMAGE10$ (x $IMAGE31$- t) $IMAGE10$/2! + … ) ∙ F(t) dt =
= K(x $IMAGE61$- x $IMAGE31$) ∙ (E $IMAGE67$F(t) dt + A∙ $IMAGE67$(x $IMAGE31$- t) ∙ F(t) dt + A $IMAGE10$/2! ∙ $IMAGE67$(x $IMAGE31$- t) $IMAGE10$ ∙ F(t) dt + … ) .
Эта формула справедлива для случая системы дифференциальных уравнений с постоянной матрицей коэффициентов A=const.
Для случая переменных коэффициентов A=A(x) можно использовать прием разделения участка (x $IMAGE61$- x $IMAGE31$) интервала интегрирования на малые подучастки, где на подучастках коэффициенты можно считать постоянными A(x $IMAGE122$)=const и тогда вектор частного решения неоднородной системы дифференциальных уравнений Y*(x $IMAGE61$←x $IMAGE31$) будет на участке складываться из соответствующих векторов подучастков, на которых матрицы Коши приближенно вычисляются при помощи формул с постоянными матрицами в экспонентах.
4 Метод «переноса краевых условий» в произвольную точку интервала интегрирования
Метод обсчитан на компьютерах. По нему уже сделано 3 кандидатских физ-мат диссертации.
Метод подходит для любых краевых задач. А для «жестких» краевых задач показано, что метод считает быстрее, чем метод С.К.Годунова до 2-х порядков (в 100 раз), а для некоторых «жестких» краевых задач не требует ортонормирования вовсе. Смотри:
Численный метод переноса краевых условий для жестких дифференциальных уравнений строительной механики
Журнал "ММ", Том: 14 (2002), Номер: 9, 3 стр. 1409-003r.pdf
Полное решение системы дифференциальных уравнений имеет вид
Y(x) = K(x←x ) ∙ Y(x ) + Y*(x←x ) .
Или можно записать:
Y(0) = K(0←x $IMAGE38$) ∙ Y(x $IMAGE38$) + Y*(0←x $IMAGE38$) .
Подставляем это выражение для Y(0) в краевые условия левого края и получаем:
U∙Y(0) = u,
U∙[ K(0←x $IMAGE38$) ∙ Y(x $IMAGE38$) + Y*(0←x $IMAGE38$) ] = u,
[ U∙ K(0←x $IMAGE38$) ] ∙ Y(x $IMAGE38$) = u - U∙Y*(0←x $IMAGE38$) .
Или получаем краевые условия, перенесенные в точку x $IMAGE38$:
U $IMAGE38$∙ Y(x $IMAGE38$) = u $IMAGE38$ ,
где U $IMAGE38$= [ U∙ K(0←x $IMAGE38$) ] и u $IMAGE38$ = u - U∙Y*(0←x $IMAGE38$) .
Далее запишем аналогично
Y(x $IMAGE38$) = K(x $IMAGE38$←x $IMAGE37$) ∙ Y(x $IMAGE37$) + Y*(x $IMAGE38$←x $IMAGE37$)
И подставим это выражение для Y(x $IMAGE38$) в перенесенные краевые условия точки x $IMAGE38$
U $IMAGE38$∙ Y(x $IMAGE38$) = u $IMAGE38$,
U $IMAGE38$∙ [ K(x $IMAGE38$←x $IMAGE37$) ∙ Y(x $IMAGE37$) + Y*(x $IMAGE38$←x $IMAGE37$) ] = u $IMAGE38$ ,
[ U $IMAGE38$∙ K(x $IMAGE38$←x $IMAGE37$) ] ∙ Y(x $IMAGE37$) = u $IMAGE38$ - U $IMAGE38$∙ Y*(x $IMAGE38$←x $IMAGE37$) ,
Или получаем краевые условия, перенесенные в точку x $IMAGE37$:
U $IMAGE37$∙ Y(x $IMAGE37$) = u $IMAGE37$ ,
где U $IMAGE37$= [ U $IMAGE38$∙ K(x $IMAGE38$←x $IMAGE37$) ] и u $IMAGE37$ = u $IMAGE38$ - U $IMAGE38$∙ Y*(x $IMAGE38$←x $IMAGE37$) .
И так в точку x $IMAGE184$ переносим матричное краевое условие с левого края и таким же образом переносим матричное краевое условие с правого края и получаем:
U $IMAGE184$∙ Y(x $IMAGE184$) = u $IMAGE184$ ,
V $IMAGE184$∙ Y(x $IMAGE184$) = v $IMAGE184$ .
Из этих двух матричных уравнений с прямоугольными горизонтальными матрицами коэффициентов очевидно получаем одну систему линейных алгебраических уравнений с квадратной матрицей коэффициентов:
$IMAGE191$ ∙ Y(x $IMAGE184$) = $IMAGE193$.
А в случае «жестких» дифференциальных уравнений предлагается применять построчное ортонормирование матричных краевых условий в процессе их переноса в рассматриваемую точку. Для этого формулы ортонормирования систем линейных алгебраических уравнений можно взять в [Березин, Жидков].
То есть, получив
U $IMAGE38$∙ Y(x $IMAGE38$) = u $IMAGE38$,
применяем к этой группе линейных алгебраических уравнений построчное ортонормирование и получаем эквивалентное матричное краевое условие:
U $IMAGE197$∙ Y(x $IMAGE38$) = u $IMAGE197$.
И теперь уже в это проортонормированное построчно уравнение подставляем
Y(x $IMAGE38$) = K(x $IMAGE38$←x $IMAGE37$) ∙ Y(x $IMAGE37$) + Y*(x $IMAGE38$←x $IMAGE37$) .
И получаем
U $IMAGE197$∙ [ K(x $IMAGE38$←x $IMAGE37$) ∙ Y(x $IMAGE37$) + Y*(x $IMAGE38$←x $IMAGE37$) ] = u $IMAGE197$ ,
[ U $IMAGE197$∙ K(x $IMAGE38$←x $IMAGE37$) ] ∙ Y(x $IMAGE37$) = u $IMAGE197$ - U $IMAGE197$∙ Y*(x $IMAGE38$←x $IMAGE37$) ,
Или получаем краевые условия, перенесенные в точку x $IMAGE37$:
U $IMAGE37$∙ Y(x $IMAGE37$) = u $IMAGE37$ ,
где U $IMAGE37$= [ U $IMAGE197$∙ K(x $IMAGE38$←x $IMAGE37$) ] и u $IMAGE37$ = u $IMAGE197$ - U $IMAGE197$∙ Y*(x $IMAGE38$←x $IMAGE37$) .
Теперь уже к этой группе линейных алгебраических уравнений применяем построчное ортонормирование и получаем эквивалентное матричное краевое условие:
U $IMAGE234$∙ Y(x $IMAGE37$) = u $IMAGE234$.
И так далее.
И аналогично поступаем с промежуточными матричными краевыми условиями, переносимыми с правого края в рассматриваемую точку.
В итоге получаем систему линейных алгебраических уравнений с квадратной матрицей коэффициентов, состоящую из двух независимо друг от друга поэтапно проортонормированных матричных краевых условий, которая решается любым известным методом для получения решения Y(x $IMAGE184$) в рассматриваемой точке x $IMAGE184$:
$IMAGE239$ ∙ Y(x $IMAGE184$) = $IMAGE241$.
5 Второй вариант метода «переноса краевых условий» в произвольную точку интервала интегрирования
Этот вариант метода еще не обсчитан на компьютерах.
Предложено выполнять интегрирование по формулам теории матриц [Гантмахер] сразу от некоторой внутренней точки интервала интегрирования к краям:
Y(0) = K(0←x) ∙ Y(x) + Y*(0←x) ,
Y(1) = K(1←x) ∙ Y(x) + Y*(1←x) .
Подставим эти формулы в краевые условия и получим:
U∙Y(0) = u,
U∙[ K(0←x) ∙ Y(x) + Y*(0←x) ] = u,
[ U∙ K(0←x) ] ∙ Y(x) = u - U∙Y*(0←x) .
и
V∙Y(1) = v,
V∙[ K(1←x) ∙ Y(x) + Y*(1←x) ] = v,
[ V∙ K(1←x) ] ∙ Y(x) = v - V∙Y*(1←x) .
То есть получаем два матричных уравнения краевых условий, перенесенные в рассматриваемую точку x:
[ U∙ K(0←x) ] ∙ Y(x) = u - U∙Y*(0←x) ,
[ V∙ K(1←x) ] ∙ Y(x) = v - V∙Y*(1←x) .
Эти уравнения аналогично объединяются в одну систему линейных алгебраических уравнений с квадратной матрицей коэффициентов для нахождения решения Y(x) в любой рассматриваемой точке x:
$IMAGE242$ ∙ Y(x) = $IMAGE243$.
В случае «жестких» дифференциальных уравнений предлагается следующий алгоритм.
Используем свойство перемножаемости матриц Коши:
K(x $IMAGE31$←x) = K(x $IMAGE43$←x $IMAGE34$) ∙ K(x $IMAGE34$←x $IMAGE36$) ∙ … ∙ K(x $IMAGE