XIII Власова Б.А., Зарубин B.C., Кувыркин Г.Н. Приближенные методы математической физики (1081417), страница 63
Текст из файла (страница 63)
Определим в этом пространстве функцию е~ = ехр(В), В Е Мн+~(К): 456 8. ОДНОМЕРНЫЕ КРАЕЕЫЕ ЭАДА ЧИ д . ехр((1+ Ь1)А) — ехр(1А) — ехр(1А) = 1пп — А ехр(1А). й ос~о Ь1 Из последнего равенства следует, что общее решение системы линейных однородных ОДУ вЂ” = -А„Т(1), определяемая Т() дй матрицей А, порядка Ю+ 1, имеет вид Т(1) = ехр( — А.1)Т, где Т вЂ” произвольный вектор из К +1.
Применяя метод вариаИ+1 аии постоянных, получаем решение системы (8.56) в матричной форме т(о= «р( — А„о2 ~.Я '1 «р(-А.(Й вЂ” ое )и, ~857) о где А. = Я 'А, причем экспоненциальную функцию от квадрат- ной матрицы А„1 можно представить сходящимся рядом (8.58) В частном случае я„= я=сопя1 > О, п = О, М, имеем А, = А/я, т.е. все собственные значения Л" = Л„,/л, т = 1, У+1, матрицы А. положительные и простые. Поэтому ехр( — А,1) = СВС, (8.59) где  — диагональная матрица порядка %+1 с элементами Л с~ ехр( — — ~, т= 1, У+1. При Х» 1 практическое использование решения непосредственно в виде (8.57) затруднено.
Рассмотрим приближенные рой матрице ехр(В) б Мм+1(К) (1Х1. Нетрудно доказать, что для любых 1, я б К и А б Мм+1(К) справедливо ехр((1+ я)А) = = ехр(1А) ехр(лА) и 8.3. Неетацнонарная задача тенлонронодное он 457 способы нахождения вектора Т(1). Для этого введем совокупность интервалов времени Мы к е М. Тогда в пределах к-го интервала Ь1ь = 1ь — Фь 1 можно представить (8.57) в виде где Тя, = Т(1я 1) -- вектор узловых значений температуры в момент времени 1ь ~.
В частности, прн 1 = 1ь для конца к-го интервалаотсюда получим Тя = Т(1ь) = ехр( — А Ыь)Тя 1+ + ехр( — А,(Ыь — т))~(1ь 1+т) йт. (8.60) о Для малых значений Ь1ь в (8.58) под знаком суммы можно ограничиться лишь первым слагаемым и вместо (8.60) записать Т„=Т,, -А„Т,,Л1„+ У(1,, +т)1то Ые — А, ~(1ь 1+ т)(Ы~ — т) Нт. (8.61) о При известной зависимости координат вектора у(1) от 1 интегралы в (8.61) нетрудно вычислить прн помощи кнадрагпдрньи: Фор.ндл. Требование малости интервала Ь1ь количественно можно охарактеризовать при помощи какой-либо нормы матрицы А., например ее спектральной нормы ((А.Й.
Тогда справедливость (8.61) связана с предположением 0А.'0Ь4 « 1. Но алгоритм вычислений по (8.61) будет устойчивым по отношению к накопле пню погрешностей при выполнении более слабого ограничения й! — А.Ь1ь)( < 1, приводящего к неравенству ()А.!)Ьгь < 2. 458 в. ОДНОМЕРНЫЕ КРАЕВЫЕ ЗАДАЧИ Спектральная норма квадратной матрицы является кольцевой ~1Ч), т.е. для спектральной нормы произведения квадратных матриц В и С справедливо неравенство !)ВС!! < )!В)!!!С)!.
Так как матрица о' диагональная с элементами — > 0 на диа— ! 1 Я гонали, то все ее собственные значения положительны (как и для матрицы А). Поэтому спектральные нормы матриц В и А будут равны их наибольшим собственным значениям [1Н]. Тогда с учетом (8.38) можно написать (!А !! = !(В А!! < ( шиах — 1 шах Л„( 2 — ", Хи=О,Мниl и=О,М яппи где 5„,.„ = п1ах 5„, =о,ч лини = пнп л„.
ипО,Ж Таким образом, иэ условия !!А„!!Ыь < 2 вычислигнельной усвпой- чивосши рассматриваемого алгоритма следует ограничение на выбор интервала времени (8.62) 'Смэ Амосов А,А., Дубппския Ю.А., Квичвнова Н.В. В том случае, когда собственные значения каждой из матриц В ' и А сильно отличаются друг от друга, т.е. их число обусловленности велики, ограничение (8.62) заставляет выбирать Ыь слишком малым и использование (8.61) становится нерациональным.
Сиспэемы ОДУ вида (8.56) с такими матрицами принято называть жесюквлвв. В рассматриваемой задаче теплопроводпостн такая ситуация может возникнуть вследствие сильной неоднородности материала тела, т.е. резкого изменения функций с(х) и А(х), н при значительной неравномерности сетки. Аналитическое решение для жестких систем содержит обычно зкспоненциальные функции с сильно отличаюшимися друг от друга показателями, что и вызывает указанные вычислительные трудности и требует применения специальных методов".
459 8.3. Неетацнонапная эадача теплопроаоднооти Для численного решения системы (8.56) можно использовать один нз вариантов метпода Рунге — Ьуэнпэы. Эти варианты можно построить на основе представления решения этой системы интегралом в пределах интервала времени Ь4: т =Т,,АА ' /(т)А А ) — Ат)А, - ))А. )8.63) о Если интеграл в (8.63) заменить с погрешностью, пропорциональной (Ыь)з, выражением ® ) — АТь ))Ыы где ~ь = у(~ь ) ), то получим рекуррентную формулу мепюди Эйлери Т, = (К- А.~ Е„)Т„, +~-'~а ) Л~,. (8.64) Метод Эйлера является простейшим вариантом метода Рунге — Кутты.
В этом случае условием устойчивости алгоритма остается неравенство (8.62). Если для представления интеграла в (8.63) использовать квадратурную формулу тприпении, то придем к выражению Ть = Ть ) — А. (Та-) + Ть) — + 5 ®-) + Уь) — (8 65) 5~а ) ЛС„ 2 с погрешностью, пропорциональной (Ыь)~. Отсюда следует, что Ть=(1+-А.Ыь) ((1 гА.Аа~кьТк-)+зо (Уа-)+Уь)Ыа) где уь = ~(~ь). Условие устойчивости алгоритма вычислений теперь принимает вид ~) (~+ -А„Ась) (1 — — А.М.~ь )! < 1. (8.66) Можно показать, что для матрицы А с частичным диагональным преобладанием н диагональными элементами йн ) 0 зто условие выполняется для любых значений Ь1ы 460 8.
ОДНОМЕРНЫЕ КРАЕВЫЕ ЗАДАЧИ Чтобы избежать обрарцения матрицы 7+ Агарь/2 и сохранить третий порядок погрешности аппроксимации, алгоритм вычислений на Й-м интервале времени можно построить следующим образом: сначала по (8.64) вычислить вектор-прогноз Ты который затем использовать в правой части (8.65) вместо вектора Ть. Условием устойчивости такого алгоритма с прогнозом и последующей коррекцией будет неравенство (8.67) которое обычно в меньшей мере ограничивает выбор Ыь, чем условие (8.62). Другой путь построения алгоритма, обеспечивающего третий порядок погрешности аппроксимации на Й-м интервале времени, связан с представлением интеграла в (8.63) выражением (~ь 171 — АТь 171)ЬРь 1, где индекс и — 1рР2 указывает на то, что значения зависящих от времени величин соответству- 1 ют моменту ря 171 =1ь 1+-ЬРЬ в середине этого интервала.
2 Тогда по методу Эйлера найдем вектор 1 х 1 Тр,-17э = (7 — — А.варь)Тк-р+ -Я р~р-рЬСр, (8 68) а затем — искомый вектор Ть = Ть-1+ 5 (7ь 171 — АТь 171) ЬРя 1. (8.69) Устойчивость двухступенчатого алгоритма (8.68), (8.69) гарантирована при выполнении неравенства (8.67). Отметим, что в соответствии с (8.58) левые части в (8.62) и (8.67) являются нормами усеченных представлений выражения ехр( — А.Ьрь) г остаточными членами 0((ЬСь)э) и 0((Ьрь)э) соответственно. Чтобы найти вектор Ть с погрешностью порядка г = 4 необходимо в пределах Й-го интервала времени использовать 8.3. Неотациоиариая задача теплоироводиооти 461 трехступенчатый алгоритм, например такой: ЬТь('1= 5 '(~д ~ — АТь г)Ыь 1 ЬТ~ ~1= 5 '~~ь,~г — А(Тк г+ — ЬТ„' ))Ь(ь, ЬТ„~=Я '(Уу,— А(Ть 1 — ЛТь +2ЬТь ))Ь|ы Ть = Ть 1+ — (ЬТ +4ЬТь +ЬТь ). 6 Нетрудно проверить, что условием устойчивости этого алго- ритма будет неравенство (А а,(„)г (А а( )з Наконец, приведем пример четырехступенчатого алгоритма с (1) (г) порядком погрешности и = 5, в котором для ЬТ~ и ЬТ~ сохраняются прежние выражения, а далее вычисления проводят по формулам 1 ЬТь — -5 ~(~ь г7г — А(Ть г+ -ЬТ ))Ь(ь, ЬТ,"1 = я-'(у, — А(т,, + ~гТ„"1)) д(,, Ть = Ть г + -(ЬТь + 2ЬТ„~ + 2ЬТь 1+ ЬТ„~).
Ть = Ть гехр( — А.Ь(ь). (8.70) Все перечисленные алгоритмы принадлежат семейству метода Рунге — Кутты. В этом семействе можно построить алгоритмы с еще более высоким порядком погрешности и. Способы их построения удобно пояснить на примере решения системы (8.56) с нулевым вектором Д1).
В этом случае на й-м интервале времени иэ (8.60) получим и. ОДНОМЕРНЫЕ КРАЕВЫЕ ЗАДАЧИ Учитывая гвойгтво экспоиеициальиой функции ехр( — А,Ь1ь) = ехр(--А.Ь1ь~ ехр( — -А.Д<1ь), можно в г<ютветгтвии г (8.58) сначала вычиглить ехр( — — А.Ь1ь) 1+ ~~< ч А'„ ( — ! )<(ЬСя/2)< з=! 3.' г порядком погрешности г, а затем последовательно за г шагов г тем же порядком погрешиогти найти ехр( —, А,<51я), ехр( — —,А.Ыь), ..., ехр( — -А,ЬХь) и, наконец, ехр( — А.11<1ь), использовав этот результат в (8.70). Наиболее употребительными являются алгоритмы метода Рунге — Нутты при г < 5.
Их можно использовать и для решения нелинейных задач (например, в случае, когда в (8.49) параметры е, Л и 1~~ произвольным образом зависят от искомой функции — температуры). Реализация третьего подхода к решению пестациоиариой краевой задачи при помощи МИР связана г предположением о постояпгтве гкорости изменения искомой функции в пределах интервала времени в каждом узле одиомериой сетки по прострапгтвеииой координате и равносильна использованию в (8.56) аппроксимации производной по1 копечно-разиостиым гоотношением вида (8.52). Тогда из (8.56) получим разиостиую гхему Т вЂ” Т Я =(1 — <7)Ць < — АТь <)+<1(~~ — АТь) (8,71) Ь1„ с весовыми коэффициентами ! — <1 и и Е (0,1], называемую двухслойной раэностиной схемой с весома. погкольку в иее входят узловые значения температуры иа двух слоях нрострвнг<ввгнно-временной сетки.
8,3. Нестапнонарпав э«дача теплопроводиости 463 Значению и = О соответствует лвмал двухслойная разнос!танал схемо, совпадающая с (8.64). В этом случае (8.71) можно явно разрешить относительно искомого вектора ть, а каждое узловое значение Т~, п = О, М, температуры в конце й-го интервала времени д1ь вычислить независимо от остальных из выражения т„"=т„"-!'(1- — ид! 1 "" "-'+'" *+'+"и д1„ ви ви которое получается из (8.53) при аппроксимации производной сЕТ„(!) ~Й правой конечной разностью (разностью вперед"). Если в (8 49) с(х) =с=сопв1, Л(х) = А =сопв1, 1~~~(1,х) =О, то в случае одномерной равномерной сетки с шагом и разностная схема, соответствующая (8.72), принимает вид т! ! 2т«!+та ! и — ! и и+! и Ти — Та ' и и ь2 Т" « Т" ««1 -1! Т" « Т" ' «-! Т" ' «-! и-1 и.!-1 л-1 Рис.
8.3 Рис. 8.4 Сюда входят узловые значения температуры в четырех узлах пространственно-временной сетки, составляющих нлоблом этой разностоной схемы. При графическом изображении шаблона его узлы принято соединять отрезками прямых (рис. 8.3). В данном случае шаблон и разностную схему называют четырехточечными.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
