Введение в вычислительную физику - Федоренко Р.П.
ISBN 5-7417-0002-0
Скачать (прямая ссылка):
Метод граничных сеток. Перейдем к описанию метода уточнения расчетов около dG. Сначала, однако, обсудим возможности стандартных способов такого уточнения. Простое уменьшение шага, как уже отмечалось, делает расчеты слишком дорогими. В технике метода конечных элементов проблемы адаптации к сложному, нерегулярному виду границы области решают, используя конечные элементы с носителями неправильной формы. Это позволяет аппроксимировать
496
ПРИБЛИЖЕННЫЕ МЕТОДЫ ВЫЧИСЛИТЕЛЬНОЙ ФИЗИКИ
[Ч. II
границу не ступенчатой кривой, а, например, набором малых хорд. Сложный характер решения вблизи границы (10) учитывается тем, что вместо стандартных базисных функций в ячейках сетки, примыкающих к границе, вводятся функции, уже содержащие особенность требуемого типа (корневую, в данном случае).
Однако, как нетрудно понять, применение таких методов разрушает «теплицеву» форму матрицы А, и приходится работать с общей четырехиндексной матрицей A1-t J1, причем эта матрица — индивидуальная для каждой области. К таким же последствиям приводит и прием «регуляризации». Можно искать решение в форме и(х, у) = ^Ф(х, у) v(x, у), где Ф(х, у) — известная функция, положительная в G и без касания обращающаяся в нуль на 6G (построение такой функции в общем случае не так-то просто!), v(x, у) — подлежащая расчету функция, уже не содержащая корневой особенности.
Изложенный ниже метод позволяет заметно уточнить расчет и около границ, достаточно надежно определить N{\) и сделать первые шаги в расчете роста трещин без существенного увеличения объема вычислений. Начнем с описания области G. Она задается своим контуром dG, послед-Рис. 56 ний — набором вершин (Xi, Yi} (і =1,2,..., /).
Каждые две соседние вершины контура (і-я и (i + 1 )-я) соединяются отрезком прямой (который называется г'-м ребром). Таким образом (так как I-я вершина совпадает с первой) мы получаем замкнутую кривую («полигон»). При достаточно большом I (около 100 в расчетах) эта кривая хорошо аппроксимирует те, в конце концов не такие уж сложные контуры, которые встречаются в прикладных расчетах. Первый этап расчета проводится так, как это было описано выше, на сетке с некоторым относительно грубым шагом А. Он дает нам «грубое решение» Ukt т, которое по интерполяционной формуле (6) восполняется до непрерывной в Gh функции U(x, у).
Следующий этап расчета — уточнение решения около контура. Он состоит из I отдельных задач. Около каждого г-го ребра строится своя локальная малая область Coi, покрытая сеткой с шагом Hi. Сетка строится так, что г-е ребро проходит по линии сетки (рис. 56), а область Oii покрывает некоторую окрестность г'-го ребра. После этого решается задача (1) в предположении, что в GNcoi решение уже известно (это «грубое» решение U). Такая локальная задача имеет форму
= у), (х, у) є COi, (11)
ш
/,(*, У) = fix, У) + ? д J J 7 U(*'> У') dx' dy'. (12)
G\a>,
§30]
ПСЕВДОДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ
497
Задача (11) решается стандартным образом. Хотя число таких уточняющих задач велико (около 100), число узлов в Oil относительно мало (около 100 при Ai « A/З). Объем вычислений так сильно зависит от числа узлов, что решение всех вспомогательных задач требует примерно того же машинного времени, что и получение грубого решения. Следует подчеркнуть, что наиболее трудоемкий элемент вышеизложенной методики — не решение системы (11), а ее формирование, т.е. вычисление /, по формуле (12). Именно этому элементу надо уделить особое внимание.
Функции /г вычислялись так. Сетка с шагом Ai продолжалась на несколько шагов за пределы шг. Вдоль границы со; выделялся некоторый «пояс», покрытый сеткой с малым шагом Ai. Интеграл по поясу от грубого решения вычислялся достаточно аккуратно с помощью матрицы AhI. В этом поясе значение U интерполировалось в узлы сетки по формуле (6). Такая аккуратная процедура связана с тем, что при вычислении этой части интеграла нужно учитывать сингулярность ядра интегрального преобразования (I). В остальной части GNcoi ситуация проще. Ядро уже гладкое, оно быстро убывает
(как 1 /г3). Функция U(x, у) тоже достаточно гладкая. Поэтому соответствующая часть интеграла в (12) вычисляется на более грубой сетке с шагом (2 -н 3) А. Разумеется, точность подобной методики нуждается в тщательном контроле.
Проверка методики. Описанный выше метод расчета трещин проходит стадию становления; опыт его применения пока ограничен. В такой ситуации естественно встает вопрос о доверии к полученным результатам, о контроле самой методики. В вычислительной физике эта проблема всегда возникает при освоении нового класса задач. И решается она специфическими, не очень строгими методами. Математически строгие оценки либо отсутствуют, либо настолько грубы, что реального представления о точности расчетов не дают. Ниже мы опишем и обсудим те средства контроля, которые использовались в этой задаче. По своему характеру они типичны в вопросах подобного рода.
Основное средство контроля — сопоставление расчетов с некоторыми известными точными решениями (например, известны решения (1) в эллипсе при линейной зависимости / от х, у). Такие сравнения проводились, их результат был признан положительным. Что это означает, мы обсуждать не будем, отсылая читателя к специальным подробным публикациям. С методической точки зрения здесь все ясно. То же самое относится и к сравнению с некоторыми приближенными аналитическими решениями (например, для трещин в форме вытянутого прямоугольника известны асимптотики решения вблизи середины длинных сторон границы).
