Введение в вычислительную физику - Федоренко Р.П.
ISBN 5-7417-0002-0
Скачать (прямая ссылка):
514
ПРИБЛИЖЕННЫЕ МЕТОДЫ ВЫЧИСЛИТЕЛЬНОЙ ФИЗИКИ [Ч. II
Итак, алгоритм построен. Он основан не на точной теории, а на соображениях, привлекаемых из близких ситуаций (в которых эти соображения являются результатом достаточно аккуратной теории). Такой способ действия характерен для вычислительной физики. В ней редко встречаются чистые, укладывающиеся в уже готовую теорию задачи. Специалист по вычислительной физике обычно начинает построение алгоритма «по аналогии» с тем, что он уже знает, и начинает именно с практического использования своего алгоритма, а не с развития соответствующей ему теории. Это понятно: ведь алгоритм может оказаться неудачным и стоит ли тогда строить теорию? Если же он оказался удачным, наступает время решать прикладные задачи, а теория может и подождать.
. В данном случае, конечно, в первую очередь хочется понять, удачен ли алгоритм или надо в нем что-то менять. Естественным средством проверки алгоритма является решение задачи, имеющей точное решение и содержащей характерные для данного случая трудности. Такая задача может быть построена. В круге единичного радиуса в центре помещена скважина малого радиуса. Коэффициент диффузии берется в виде (1 — г)312. Поскольку задача цилинд-рически-симметрична (ее решение зависит только от г), уравнение диффузии становится обыкновенным уравнением:
1=1? [r^1 “ r)3/27F] =Z = cOnst, г Є [0,1].
Оно элементарно интегрируется:
,_____ / + 2 С. і JT~~f
„(,)_/VT=T +-J=I+C1 In +C1.
где C1, C2 — произвольные постоянные.
Ограниченное решение находим при C1 = —//2. Постоянная C2 не существенна, положим C2 = 0. Это решение имеет корневую особенность на внешней границе (г = 1) и логарифмическую — в центре. Преобразованием подобия получим решение уравнения (2) в области р =? г(х, у) < R, ще р = 0.1, R = 250. Постоянную / легко подобрать так, чтобы выполнялось внутреннее краевое условие: р = 1.2 на границе скважины. Задача решалась с помощью метода конечных суперэлементов на квадратной сетке с шагом H — 17. Число счетных узлов было около 650. Заметим, что граница исходной области (окружность г = R) аппроксимируется контуром (Zi, Yi] достаточно аккуратно, хбтя область определения счетных величин Pk т (множество счетных узлов) аппроксимирует круг г < R очень грубо.
Вычислительный эксперимент имел целью выяснить два обстоятельства: как сходятся итерации и какова точность разностного решения? Сходимость итерационного процесса (использовался второй способ) иллюстрирует табл. 24, в которой представлены: v — число итераций, є — невязка (максимум по (к, т) модуля левой части (16)),
§31]
МЕТОД КОНЕЧНЫХ СУПЕРЭЛЕМЕНТОВ
515
Таблица 24
О) V є р' р"
1 200 6.5-IO-4 0.568 -0.081
1 400 I О <Г5 0.626 0.018
1 700 1.0- IO-4 0.661 0.079
1.5 400 1.8- IO-5 0.675 0.103
1.7 200 2.8-IO-5 0.675 0.101
1.75 400 4.6-IO-8 0.678 0.107
1.8 400 1.6-IO-9 0.678 0.107
со — значение параметра релаксации в (18). Итерации начинались с Pk т = 0, значения ск т « 4 -S- 5. Из таблицы видно, что при со = 1 невязка быстро достигает малых значений, затем сходимость становится очень медленной: за первые 200 итераций невязка уменьшается почти в IO4 раз, в дальнейшем за 200 итераций она уменьшается примерно вдвое.
Как расценить этот результат? С одной стороны, погрешность в правой части (порядка
0.01 %), кажется, не требует существенного улучшения, с другой стороны, медленная сходимость внушает какую-то тревогу.
Обычно она свидетельствует о том, что оператор (здесь div и3 grad) имеет собственное значение, близкое к нулю. В § 14 специально отмечалось, что связь между невязкой и погрешностью определяется минимальным собственным значением: погрешность есть величина порядка невязки, деленной на | Xmin |.
В рассматриваемой задаче есть основания подозревать наличие очень малого собственного значения. В самом деле, легко угадать функцию, которая является «почти собственной» с собственным числом Х = 0. Это есть функция р(х, у) = 1. Она удовлетворяет «почти
всем уравнениям» div (и3 gradp) = 0. He выполнено только однородное краевое условие на скважине: р = 0 на Sg. Конечно,
X = 0 не является собственным значением, но приведенное выше рассуждение служит основанием ожидать близкого к нулю собственного значения, тем более близкого, чем меньше радиус скважины.
Действительно, результаты, приведенные в табл. 24, показывают, что дальнейшее уменьшение невязки (казалось бы, излишнее) сопровождается заметным изменением важных величин р', р". Это значения рк т, взятые на луче, выходящем из центра скважины по диагонали сетки: р' — первое значение р на луче (на расстоянии Я/VT от центра скважины), р" — последнее, почти граничное значение р на луче. Видно, что уточнение в процессе итераций решения системы разностных уравнений до значений є « 10_6 имело смысл.
Таблица 24 дает представление и о точности самой разностной схемы (при шаге H= 17). Значения точного решения в соответст-
Рис. 66
516
ПРИБЛИЖЕННЫЕ МЕТОДЫ ВЫЧИСЛИТЕЛЬНОЙ ФИЗИКИ
[Ч. II
вующих точках суть 0.680 и 0.124. На рис. 66 представлено точное решение в сечении ПО ЛИНИИ х = у, проходящей через центр скважины. Кстати, значение р на скважине есть 1.2. Приближенное решение отмечено кружками.
