Ортега Дж., Пул У. Введение в численные методы решения дифференциальных уравнений. Под ред. А.А.Абрамова (1986) (1095855), страница 57
Текст из файла (страница 57)
Если вычитание из А матрицы (Ьт/2) В не приведет к существенному осложнению реализации У. 1У-разложения, то затраты при счете по методу (7.4.13) не будут намного большими, чем в методе Эйлера. Однако, как мы видели в гл. 2, метод трапеций имеет второй порядок точности и позволяет использовать для жестких уравнений более крупный шаг по времени. Следовательно, этот метод больше подходит для данной задачи. Тот же самый подход можно применить и к гиперболическим уравнениям. Мы сейчас покажем это на примере волнового уравнения и„=си»», 0<х<1, у>0, (7.4.14) с начальными и граничными условиями и(0,х)=У(х), и,(О,х)=я(х), и(т, 0)=0, и(у,1)=0.
(7.4.15) Гиперболические уравнения, такие, как (7 4.14), обычно сводят к системе уравнений, в которую входят только первые производные по времени. Это совершенно аналогично ситуации с обыкновенными дифференциальными уравнениями. Напомним, что в гл. 2 мы сводили уравнения более высокого порядка к системзм обыкновенных дифференциальных уравнений первого порядка.
В случае (7.4.14) мы осуществим это сведение, введя функцию и(у,х) такую, что и, = аи», и, = аи„, (7.4.16) где а =ус. Если функции и и и являются решением системы(7.4.16) и обладают достаточной гладкостью, то, дифференцируя первое из уравнений (7.4.16) по т, а второе — по х, мы получаем и~с =ад„г = акр» = а и»» = си»», т.е. и — решение (7.4.14). При этом начальные и граничные условия (7.4.15) принимают вид и а й(х, у) = Х ау(у)ру(х), и(х„г) = Х Д(т)$;(х). (7.4.19) а=т Если мы потребуем, чтобы эти приближенные решения удовлетворяли уравнениям (7.4.16) в точках х „..., х„, то придем к соотношениям и а Х ау(г)ру(ху)=а Х Д(т) Ф,".(ху), у=],...,п, (7.4.20) Х Р,'(г)Фг(ху) =а Х а;Яр,'(ху), у =1,...,и, у= 1 !=1 ]» и(О,х) =у(х), п(О,х) = — у" у(т)сй, а 0 и(т, 0)=0, и(т, 1)=0. (7.4Л7) Пусть теперь р1 (х),... „р „(х) и Ф, (х),..., Д„(х) — два набора базисных функций, удовлетворяющих условиям ~ру(0) =~ру(1) = ~;(0) = 4(у(1) — О, г =!,..., л.
(7.4.18) ~ л л Мы будем искать приближенные решения и и п системы (7.4,16) в виде представляющим собой связанную систему обыкновенных дифференциаль- ных уравнений относительно неизвестных функций а1...а„и 01,... ..., Р„. Как и раньше, из (7.4.17) получаем начальные условия и Х а~(0) у;(х;) = 1'(х;), 1=1 и х~ Х Р,.(0)~,(х;) = — ( ()Л, ~=1 а е Мы здесь опять предполагаем, что матрицы (~рр (х )) и ф~(х;)) являются невырожденными.
Таким образом, полудискретный метод для волнового уравнения (7.4.14) полностью аналогичен соответствующему методу для уравнения теплопроводности с единственным исключением, что система обыкновенных дифференциальных уравнений теперь содержит вдвое боль- ше неизвестных функций. В двух приведенных примерах мы при дискретизации по пространствен- ной переменной основывались на принципе коллокации, хотя могли бы воспользоваться также и рассмотренным в гл, 5 методм Галеркина. Для дискретизации можно использовать и конечно-разностную аппрок- симацию, которую сейчас кратко обсудим на примере уравнения тепло- проводности (7.4.1) .
Если и — точное решение (7.4.1), то в узловых точ- ках х,,..., х„должны выполняться приближенные соотношения 1=1,...,и, с и~(1, х;) [и(г, х;+1) — 2и(г,х~)+и(г, х; 1)]. (7.4.21) Это ведет к следующей процедуре. Мы хотим найти и функций и, (г),... ..., с„(г) таких, что и;(г) и(г, х;), ю'= 1,..., и. Приближенные соотношения (7.4,21) наводят на мысль, что эти функции можно попытаться найти как решение системы обыкновенных дифференциальных уравнений с ир.(г) = — ~ин 1 (г) — 2и~(г) + с~ 1 (г)3, (7.4.22) й~ в которой в соответствии с граничными условиями (7.4.2) функции ие и и„+, считаются тождественно равными нулю. Кроме того, из начального условия (7.4.2) мы имеем и;(О) = У'(х;), с = 1,..., и. (7.4.23) Система (7.4.22) удобно записывается в матричной форме как т (г) = А т(г), (7,4,24) (Ьх) где теперь А — просто трехдиагональная матрица (7,2.14).
Если мы применим к этой системе метод Эйлера, то получим сЬ| т"'+'=т — — Ате, т =О, 1,... (Ьх)2 В покоординатной форме мы приходим к формулам сЬг т+1 т+ ( т 2 тт + т 1 1 (ЬХ)2 Н гп=0,1,..., представляющим собой просто явный метод (7.25). Аналогично этому неявный метод (7.3.2) получается при применении к (7.4.24) обратного метода Эйлера, а метод Кранка — Николсона (7.3.12) возникает в результате применения к (7.4.24) метода трапеций (2.5.38) . Описанную процедуру, приводящую к системе обыкновеных дифференциальных уравнений (7.4.24), иногда называют методом прямых, и потенциально этот метод имеет очень широкую область применения.
~г~ополнигельные замечания и ссылки 7.4 Вопросы, связанные с использованием полудискретных методов, являются предметом активного исследования, и пока еще нельзя окончательно сказать, сколь велика эффективность этих методов по сравнению с другими. Дальнейшее обсуждение этих вопросов можно найти в книге [47] . Наиболее универсальное математическое обеспечение численного решения дифференциальных уравнений в частных производных, включаюших зависимость от времени, основывается на методе прямых (см. [42]).
Преимущество этого подхода заключается в возможности использовать отличное математическое обеспечение, разработанное для решения обыкновенных дифференциальных уравнений [96]. Описание основанных на методе прямых различных пакетов математического обеспечения для решения дифФеренциальных уравнений в частных производных общего вида имеется в работах [43, 45, 72, 104]. УПРАЖНЕНИЯ 7.4 7.4.1. а) Выпишите явно систему уравнений (7.4.7) для базисных функций еу,(х) = = з1п(х их), й = 1... л, считая, что узловые точки х„...,х„нахоцятся наодинаковом расстоянии друг от друга. Составьте программу, реализующую метод Эйлера(7.4.12) при и = 10 и Ь1 = 0,1 Просчитайте по этой программе 20 шагов по времени, задавая различные начальные условия; б) проделайте то же самое, используя в качестве ру, квадратичные сплайны, рассмотренные в разделе 5.2.
7 42. Для системы (7.4.8) выпишите в явном виде формулы методов Рунге — Кутта и Адамса — Башфорта второго и четвертого порядков. 7.4.3. Составьте программу, реализующую метод трапеций (7.4.13) для задач нз упражнения 7.4.1. 7.4.4- Проделайте снова упрамснение 7.4.1 для уравнений (7,4.20), возникающих из волного уравнения.
(.читайте, что в пункте а) все у1 и чч являются тригонометрическими функциями, а в пункте б) — квацратичными сллайнами. Глава 8 ПРОКЛЯТИЕ РАЗМЕРНОСТИ 8.1. Задачи с двумя и тремя пространственными переменными В предыдущей главе мы рассматривали уравнения в частных произвоцных с двумя независимыми переменными: временем и одной пространственной переменной.
Поскольку физические процессы происходят в трехмерном мире, математические модели лишь с одной пространственной переменной связаны с существенным упрощением реальной физической ситуации. Правда, зтн модели часто оказываются достаточными для таких физических явлений,' которые обладают различными снмметрнямн нли в которых изменения по двум из трех измерений происходят настолько медленно, по ими можно пренебречь. Однако крупномасштабные задачи научного программирования сейчас все в большей степени включают более детальный анализ проблем, в которых интерес представляют все три пространственных направления или по крайней мере два из них. В этой главе мы будем иметь дело с задачами, которые описываются уравнениями с несколькими пространственными переменными, однако для простоты изложения мы ограничим любые детальные рассмотрения двумерным случаем. В предьщущеи главе мы вывели уравнение теплопроводности (8.1.1) и, =сига как математическую модель изменения температуры в длинном тонком стержне.
Давайте теперь рассмотрим ту же самую задачу для трехмерного куба, изображенного на рис. 8.1. Вывод уравнения (8.1.1) непосредственно распространяется на случай трех измерений, и при этом в уравнении появляются частные производные по всем трем переменным х,унт. Таким образом, обобщение (8.1.1) на случай трех пространственных переменных имеет вид иг =с(ихх +иуу +игг) (8.1.2) где по-прежнему константа с есть отношение коэффициента теплопроводности к произведению удельной теплоемкости на плотность.
(Предполагается, что эти величины являются постоянными, т.е. не зависят от времени и пространственных переменных.) Уравнение (8.1.2) моделирует температуру и как функцию времени и точек внутри тела. Для завершенности модели мы, как обычно, должны задать начальные и граничные условия. Задание граничных условий проще всего продемонстрировать на соответствующем двумерном случае: и, = с(и„„+иуу). (8.1.3) Уравнение (8.1.3) можно рассматривать как математическую модель 1б. Дж.
Ортега 241 изменения температуры в тонкой плоской пластинке. Мы будем считать, что эта пластинка представляет собой единичный квадрат, показанный на рис. 8.2. Простейшим видом граничных условий является задание температуры на всех четырех сторонах пластинки: и (~, х, у) = К(х, у) на сторонах пластинки, (8.!.4) где 8 — заданная функция. В качестве другой возможности можно, например, предположить, что одна из сторон, скажем х = О, абсолютно изолиро- Рис. 8.1. Трехмерный куб Рис.
8.2. Плоская тонкая пластинка (8.1.6) и (О, х, у) = Дх, у). Из физических соображений ясно, что, каково бы ни было начальное распределение температуры (8.1.6), при заданных граничных условиях вида (8.1.4) и (илн) (8.'1.5) со временем должно установиться некоторое окончательное стационарное распределение, которое определяется только этими граничными условиями. Во многих случаях основной интерес представляет определение именно этого стационарного состояния. Так как стационарное решение не зависит от времени, то оно должно удовлетворять уравнению(8.1.3), где и, = О, т.е. должно удовлетворять уравнению Лапласа и„„+и„=О„ (8.1.7) которое, как упоминалось в предыдущей главе, является примером эллиптического уравнения.
Если бы нас в задаче о распределении температуры интересовало только стационарное решение, то мы, в принципе, могли бы действовать двумя способами — либо решать уравнение (8.1.3) относитель- вана. Тогда поток тепла через эту сторону равен нулю, что приводит к граничному условию и„(1, О,у) =О, О~ у<1, (8.1.5) которое может сочетаться с заданием условий (8.1.4) на других сторонах.
Граничные условия вида (8.1.5) обычно называют условиями Неймана, а вида (8.1.4) — условиями Дирихле. Возможны также различные комбинации таких условий, включая задание отличного от нуля потока тепла через границу. В трехмерном случае граничные условия задаются совершенно аналогично. Кроме граничных условий, мы должны также указать распределение температуры в некоторый начальный момент времени, которому у нас соответствует г = О. Для (8.1.3) начальное условие имеет вид ятт1 Е= (х — х,,у — у,, г — а,), (8.1.9) где а — гравитационная постоянная„а г — евклидово расстояние между точками расположения масс*). Как легко проверить (упражнение 8.1.1), вектор Г является градиентом потенциала — атт1 р(х,у, г) = г (8.1.10) т.е.