Як відбувається побудова розрахункової сітки для лінійних статичних завдань, Блог COMSOL

Як відбувається побудова розрахункової сітки для лінійних статичних завдань

Walter Frei 22/10/2013

У цій статті ми розповімо про особливості побудови сітки для лінійних статичних завдань, які вирішуються методом кінцевих елементів. Це перша публікація із ряду про методи побудови сітки. Серія статей написана, щоб дати рекомендації, які допоможуть упевнено будувати розрахункові сітки для моделей кінцевих елементів.

Побудова сітки кінцевих елементів

Розрахункова сітка кінцевих елементів необхідна двох цілей. Перша - розбиття змодельованих у САПР геометрій на менші частини чи елементи. За ними можна записати систему рівнянь, що описує рішення головного рівняння. Ще сітка використовується для відображення сфери рішень фізичних завдань. Існує похибка, пов'язана з дискретизацією геометрії та дискретизацією рішення. Розглянемо кожен випадок окремо.

Дискретизація геометрії

Розглянемо два простих геометричних тіла: куб та циліндричну оболонку.

Існує чотири різні типи елементів, які можна використовувати для побудови сіток. Це геометричні тіла - тетраедри (чотиригранники), гексаедри (шестигранники), трикутні призматичні елементи (призми) та піраміди.

Сірі кружки відображають кути або вузли елементів. Ці чотири елементи можуть бути використані у будь-якій комбінації. (Для двомірного моделювання доступні трикутні та чотирикутні елементи.) Після розгляду ілюстрації стає видно, що розрахункові сітки обох геометричних тіл можуть бути побудовані з використанням лише одного шестигранника, двох призм, трьох пірамід або п'яти чотиригранників. Як ми дізналися з раніше запису в блозі про рішення лінійнихстатичних завдань шляхом кінцевих елементів, ви завжди отримаєте рішення вже у першій ітерації методу Ньютона-Рафсона. Це правильне всім лінійних завдань кінцевих елементів незалежно від виду розрахункової сітки. Тож погляньмо на найпростішу сітку, яку ми можемо накласти на вибрані структури. Нижче ілюстрація — це графічне представлення дискретизації наших геометричних тіл на один шестигранник:

Вочевидь, що розрахункова сітка куба — це ідеальне уявлення вихідного геометричного тіла, тоді як сітка циліндричної оболонки проявляється досить погано. Фактично вона з'являється у такому вигляді тільки при виведенні зображення. Елементи завжди виводяться на екран з прямими ребрами, (це зроблено підвищення графічної продуктивності) але для дискретизації геометрії (і рішення) COMSOL зазвичай використовує елементи Лагранжа другого порядку. Отже, хоча ребра елемента завжди з'являються як прямих, всередині середовища вони представлені як:

Білі кружки відображають серединні вузли ребер елементів другого порядку. Лінії, що визначають ребра елементів, проходять через три точки. Ребра апроксимовані поліноміальним наближенням. У центрі кожної з чотиристоронніх поверхонь і в центрі об'єму шестигранних елементів Лагранжа другого порядку також є додаткові вузли (опущені для простоти). Видно, що ці вузли добре служать уявлення викривлених кордонів елементів. Більшість фізичних завдань COMSOL використовує елементи другого порядку за замовчуванням. Два винятки з цього завдання, пов'язані з перенесенням хімічних речовин і полем потоку рідини. (Поки в проблемах такого типу в основному йдеться про конвекцію, головні рівняння краще вирішувати з елементами першого порядку.) Також доступні і більш елементивисокого порядку. А ось елементи другого порядку – це відмінний компроміс між точністю та вимогами до обчислювальної продуктивності.

Ілюстрація нижче показує похибку геометричної дискретизації при побудові сітки для дуги 90° у кількості елементів першого і другого порядку:

Висновок, який можна зробити з усього вищесказаного: для зменшення похибки геометричної дискретизації рівня нижче 1% потрібно хоча б два елементи другого порядку або вісім елементів першого порядку. Фактично, два елементи другого порядку вносять похибку геометричної дискретизації менш ніж 0.1%. Густіші сітки будуть точніше представляти фігуру, але й використовують більше обчислювальних ресурсів. Це дає нам дві добрі практичні поради:

При використанні елементів першого порядку підганяйте сітку таким чином, щоб припадало хоча б 8 елементів на дугу 90° градусів
При використанні елементів другого порядку на дугу в 90° повинно бути два елементи

Тепер із цими правилами, отриманими досвідченим шляхом, ми можемо оцінити похибку, внесену під час побудови сітки. Це можна зробити з певною впевненістю навіть до безпосереднього рішення моделі. А зараз давайте звернемо нашу увагу на те, як сітка дискретизує рішення.

Дискретизація рішення

Розрахункова сітка кінцевих елементів також використовується для представлення області рішень. Рішення обчислюється за вузловими точками. При його інтерполяції використовується поліноміальний базис відновлення повного поля рішень. При розрахунку лінійних завдань методом кінцевих елементів у нас завжди вийде рішення, і не важливо, наскільки сітка груба. Вона може бути занадто точної. Щоб зрозуміти, як щільність сіткивпливає на точність рішення, поглянемо на просте завдання передачі тепла на попередні тіла:

Різні температури прикладені до протилежних граней куба та циліндричної оболонки. Термічна провідність є константою і решта поверхні ізольовані.

Рішення для випадку з кубом таке. Поле температури змінюється лінійно у тілі куба. Для обчислення правильного рішення цієї моделі буде насправді достатньо одиночного шестигранного елемента першого порядку. Звичайно, так везтиме не часто.

А тепер давайте подивимося на більш складний варіант. Ми вже бачили, що модель циліндричної оболонки буде з похибкою геометричної дискретизації через викривлені ребер, так що ми запустимо її хоча б з двома елементами другого порядку (або восьмью першого порядку) вздовж викривлених ребер. Якщо ви уважно подивіться на зображення вище, побачите, як на межах викривлені ребра елементів. Внутрішні елементи мають прямі ребра. Поки поле температури не змінюється вздовж осі циліндра, ми можемо використовувати одиночний елемент. Однак у радіальному напрямку від внутрішньої до зовнішньої поверхні нам необхідно мати достатньо елементів для дискретизації рішення. Аналітичне рішення цього випадку виглядає так: \ln(r) . Його можна порівняти з рішенням методом кінцевих елементів. Поки поліноміальні базисні функції не можуть дати ідеального опису, давайте виведемо зображення похибки у вирішенні кінцевих елементів для лінійного та квадратичного елементів:

З цього графіка стає ясно, що в міру збільшення кількості елементів моделі похибка зменшується. Це фундаментальна властивість методу кінцевих елементів: що більше елементів, то точніше рішення. Звичайно, зцим пов'язані витрати. Для вирішення великих моделей потрібно більше обчислювальних ресурсів та часу. Тепер ви помітите, що тут немає одиниць виміру по осі x. Розмір вибирається у кожному даному випадку. Коефіцієнт, при якому похибка зменшується із збільшенням густоти сітки, буде своїм для кожної моделі. Він залежить багатьох чинників. Єдиний важливий момент це те, що похибка завжди монотонно зменшуватиметься для коректно поставлених завдань.

Також ви помітите, що після певної точки похибка знову починає збільшуватися. Це відбувається коли окремі елементи сітки стають дуже маленькими, і ми стикаємося з межами числової точності. Тобто числа в нашій моделі менше, ніж можуть бути представлені комп'ютером. Ця проблема притаманна всім обчислювальним методам, як із використанням кінцевих елементів. Комп'ютери не здатні точно відображати дійсні цифри. Точка, в якій похибка починає зростати, буде приблизно такою: \sqrt> \approx 1.5 \times 10^ . Ми часто говоримо, що мінімальна досяжна похибка дорівнює 10 -6, щоб не уникати достовірної практики. Таким чином, проінтегруємо обчислену різницю між точним та розрахунковим рішенням по всій моделі:

Ми говоримо, що зазвичай у межах збільшення густоти сітки похибка \epsilon може бути приведена до такого невеликого значення, як 10-6. У будь-якому випадку на практиці ж вступні дані наших моделей часто матимуть набагато більшу невизначеність, ніж ця величина. Також пам'ятайте, що загалом ми не знаємо точного рішення. Тому ми повинні порівняти розрахункові рішення по сітках з різними розмірами та спостерігати, до яких значень сходиться результат.

Ущільнення адаптивної сітки

Хотілося б завершитицю статтю описом найкращого способу ущільнити сітку. Графіки вище показують, що похибка йде на спад, коли всі елементи моделі стають меншими. Як би там не було, в ідеальному випадку ви зробите елементи меншими лише в ділянках, де похибка висока. Середовище COMSOL реалізує це через Adaptive Mesh Refinement (ущільнення адаптивної сітки), яке спочатку створює початкову сітку. Потім середовище за ітераціями додає елементи до дільниць, де за оцінкою має бути висока похибка і після цього вирішує модель знову. Це можна робити стільки ітерацій, скільки хочеться. Ця функція працює з трикутними елементами у двовимірному просторі та чотиригранними у тривимірному. Давайте перевіримо це в контексті простого будівельного завдання з механіки - пластина з одновісним натягом і отвором, як показано на малюнку нижче. Через те, що пластина симетрична, нам треба вирішити лише одну чверть моделі.

Розраховані поля зміщення та результуючі натяги цілком однорідні на деякій відстані від отвору, але змінюються у безпосередній близькості від нього. Малюнок нижче показує початкову сітку, і навіть результати кількох ітерацій згущення разом із розрахованим полем напруги.

Зверніть увагу, як вибірково середовище COMSOL додає невеликі елементи навколо отвору. Це не повинно дивувати вас, адже нам відомо, що тут буде сильніша напруга навколо отвору. На практиці для вибору відповідної сітки рекомендується використовувати комбінацію згущення адаптивної сітки, інженерної оцінки та досвіду.