РОЗДІЛ 2
Адаптивні методи вкладених сіток
2.1. Основні положення
Метою використання адаптивного методу є:
одержання розв’язку, похибка якого не перевищує наперед заданого значення;
мінімізація машинного часу за рахунок скорочення кількості вузлів сітки та
оптимального їх розташування .
Для цього метод в залежності від поведінки розв’язку автоматично вибирає розмір
кроку по усім аргументам, тобто забезпечує пошук невідомої функції на змінній
по простору та часу сітці. На кожному кроці по часу при розв’язанні
нестаціонарних задач необхідно виконати такі дії (рис.2.1):
оцінити похибку результату. Якщо вона перевищує допустиму, результати
останнього часового кроку анулювати;
вибрати нові величини кроків по часу та просторовим координатам таким чином,
щоб при максимальній величині кроків похибка знаходилась в заданих межах.
Таким чином, основними моментами при побудові адаптивних сіток є критерій
оцінки точності та механізм визначення нового положення її вузлів при переході
на наступний часовий шар.
В якості критерію оцінки точності візьмемо величину локальної похибки для
кожного вузла, яка визначає, збільшити чи зменшити величини кроків у його
околі. Для оцінки цієї похибки одержаний на кожному кроці наближений розв’язок
можна, наприклад, порівняти з однією з наступних величин:
наближеним розв’язком, обчисленим з іншою сіткою;
наближеним розв’язком, обчисленим іншим методом;
наближеним розв’язком, екстрапольованим по попереднім часовим шарам.
Оскільки окрім порівняння фактичної похибки з допустимою необхідно з’ясувати
поведінку функції у кожному координатному напрямку, доцільно використати перший
з перелічених способів.
Таким чином, основними етапами запропонованого адаптивного методу на кожному
часовому кроці є:
обчислення розв’язку на декількох вкладених сітках;
використання одержаної інформації для оцінки похибки та визначення характеру
трансформації сітки для кожного вузла;
побудова нової сітки на основі її функцій трансформації;
апроксимація розв’язку на нову сітку.
2.2. Контроль похибки результату
Розглянемо запропонований адаптивний метод на прикладі нестаціонарної задачі
для n-вимірного рівняння теплопровідності з постійними коефіцієнтами [150]:
. (2.1)
Нехай для розв’язання диференційної задачі використовується деякий різницевий
метод з точністю апроксимації на рівномірних сітках
, (2.2)
де t - крок по часу, ha - кроки по просторовим координатам. В роботах
[133,137,148,149] показано, що така точність апроксимації зберігається і для
нерівномірних сіток при невеликій різниці між двома сусідніми кроками, причому
в якості t та ha виступають їх нормовані значення.
Як критерій контролю похибки для кожного моменту часу візьмемо локальну похибку
розв’язку , тобто норму різниці між точним та наближеним розв’язками для даного
моменту часу tk+1 при умові, що для попереднього моменту tk ці розв’язки
співпадають: . Ця величина фактично означає нову похибку, яку на кожному кроці
по часу вносить різницевий метод. Можна показати, що рiзницевi схеми, для яких
виконується рiвнiсть (2.2), мають у кожному вузлi локальну похибку
, (2.3)
де величини мiстять iнформацiю про старшi похiднi невiдомої функцiї у цьому
вузлi. Визначення у кожному вузлi дозволить оцiнити похибку результату та
вибрати нову сiтку. Отже, для моменту часу для кожного вузла необхiдно знайти
невiдому , тобто побудувати вiдповiдну систему рiвнянь.
Нехай розв'язок для моменту часу вже знайдений i треба знайти його для моменту
часу , причому перевiрити його точнiсть та вибрати нову сiтку для переходу на
наступний шар . Сiтку для шару вважаємо вже вiдомою. Якщо за допомогою
рiзницевої схеми перейти з шару на шар з базовим кроком та просторовими
кроками, що задаються вiдомою нерiвномiрною сiткою в кожному напрямку, то для
вузла з координатами виконується рiвнiсть (2.3). Якщо повторно перейти з шару k
на шар k+1 з кроком , двiчi застосувавши вiдповiдну рiзницеву схему, то
одержимо результат , для якого виконується рiвнiсть
. (2.4)
Наступнi n рiвнянь шуканої системи одержимо, якщо повторимо перехiд з k-го на
(k+1)-й шар почергово з удвiчi щiльнiшою сiткою спочатку в напрямку x1, потiм
x2 i так далi. Одержанi в кожному випадку для кожного вузла значення функцiї
позначимо через , і т.д. В результатi для кожного вузла сiтки на (k+1)-му шарi
одержимо систему лiнiйних алгебричних рiвнянь розмірністю вiдносно .
, (2.5)
Після розв'язання цієї системи за формулою (2.3) можна оцiнити для кожного
вузла породжену на даному кроцi похибку . Максимальне її значення по усім
вузлам будемо вважати локальною похибкою на (k+1)-му шарi і вона порівнюється з
допустимою похибкою . Якщо , то результати анулюються і на базi k-го шару знову
шукаємо розв'язок на ущільненій сітці. Якщо , крок приймається і формула (2.3)
дає змогу уточнити результат:
.
2.3. Визначення функцій трансформації сітки
Незалежно від того, чи буде виконуватися умова , наявна інформація про
поведінку функції дає можливість побудувати нову сітку. Вона будується для
моменту часу у випадку відмови () і для моменту , якщо крок приймається (). Для
вибору нової сiтки знайдемо допустимi коeфiцiєнти змiни кроку в околi кожного
вузла сітки на (k+1)-му шарі (), тобто функції трансформації сітки. Очевидно,
що при змiнi крокiв значення прогнозованої похибки не повинно перевищувати
eдоп:
або
. (2.6)
Оскільки кожен з доданків цього виразу означає складову похибки в якомусь
координатному напрямку, доцільно вимагати, щоб вони були рівними
. (