Алгоритмічні аспекти методу скінченних елементів

В главі І розглянутий комбінований метод для розв'язування скінченоелементної апроксимації задачі про згин пластинки, що базується на маршевому методі і методі блочної верхньої релаксації та є більш ефективним у порівнянні з іншими методами.

Розглянемо задачу згина пластики в прямокутнику

u⁽⁰⁾ = f(x);

.

Для розв'язання цієї задачі скористуємось наступною різницевою схемою В.Г. Корнєєва. скінченний елемент лаплас передобумовлювач

6/h³ [(2/h (-u_{і-1, j} ⁽⁰⁾ + 2u_{і, j} ⁽⁰⁾ - u_{і+1, j} ⁽⁰⁾) - (u_{і-1, j} ⁽¹⁾ + u_{і+1, j} ⁽¹⁾) +

2/h(-u_{і, j-1} ⁽⁰⁾ + 2u_{і, j} ⁽⁰⁾ - u_{і, j+1} ⁽⁰⁾) - u_{і, j-1} ⁽²⁾ + u_{і, j+1} ⁽²⁾)]= f(x^і).

1/h² [(6/h (u_{і-1, j} ⁽⁰⁾ - u_{і+1, j} ⁽⁰⁾) + 2(u_{і-1, j} ⁽¹⁾ + 4u_{і, j} ⁽¹⁾ + u_{і+1, j} ⁽¹⁾) +

+ (-u_{і, j-1} ⁽¹⁾ + 2u_{і, j} ⁽¹⁾ - u_{і, j+1} ⁽¹⁾)]= 0.

1/h² [(6/h (u_{і-1, j} ⁽⁰⁾ - u_{і, j+1} ⁽⁰⁾) + 2(u_{і, j-1} ⁽²⁾ + 4u_{і, j} ⁽²⁾ + u_{і, j+1} ⁽²⁾) +

+ (-u_{і-1, j} ⁽²⁾ + 2u_{і, j} ⁽²⁾ - u_{і+1, j} ⁽²⁾)]= 0.

Тут u_h⁽⁰⁾, u_h⁽¹⁾, u_h ⁽²⁾ - сіткові функції, які відповідають значенням u_h⁽⁰⁾, du⁽⁰⁾/dх ₁, du⁽⁰⁾/dх ₂, що задані в вузлах х^і прямокутної сітки, причому u^і = 0, якщо х^і , і= (і₁, і₂, і₃). Систему різницевих рівнянь запишемо у вигляді:

My = f,

,

С = В^Т,

,

,

,

.

Введемо наступне сімейство матриць:

Р_-1 = 0, Р ₀ = 1,

В Р_k = A P_k-1 - C Р_k-2, k = 1, …, n,

R_-1 = 0, Р ₀ = 1,

C R_і = A R_і-1 - B R_і-2 і = 1, …, n,

(1)

(2)

(3)

(4)

Теорема 1. Нехай у - розв'язок системи рівнянь М y = f, тоді для будь-яких j, q, k три вектори:

y_j, y_q, y_k, 0 j < q k n. y₀ = y_n+1 = 0,

задовольняють співвідношенню:

- C (R_q-j-1)^-1 y_j + Q_k^j+1 (q) y_q - B (P_k-q)^-1 y_k+1 = C (R_q-j-1)^-1 V^q_q-j + B (P_k-q)^-1 u^q+1 _k-q (5).

Нехай існують цілі позитивні числа s і p такі, що:

(p + 1) s = n + 1,

де n - блочний порядок матриці М. Тоді для побудови редукованої системи виберемо індекси j* наступним чином: j*=0, s, 2s, …, (p + 1) s. Введемо позначення:

x_l = y_sl, 0 l p + 1,

де x₀ = y₀ = 0, x_p+1, = y_n+1 = 0.

Використовуючи теорему 1:

з j = (l - 1)S, q = ls, k = (l + 1)S = 1,

одержуємо редуковану систему відносно невідомих х ₁,…х_р:

- G_l x_l-1 + Q_l x_{l -} H_l x_l+1 = F₁, l= 1, …, p. (6)

Отриману систему (6) будемо розв'язувати блочними ітераційними методами (метод Якобі або блочний метод верхньої релаксації). Розглянемо блочний ітераційний метод:

(7)

G_l = C (R_S-1)^-1, H_l = B (P_S-1)^-1,

.

Теорема 2. Ітераційний процес (7) для системи М y = f збігається.

Для оцінки алгоритмічної складності обчислень при розв'язуванні даної задачі в § 1.4 розглядається достатньо загальний оціночний функціонал, що дозволяє отримати оцінку різноманітних характеристик обчислювального процесу: числа арифметичних дій, числа макрооперацій, об'єму необхідної пам'яті тощо. При певній інтерпретації цього функціоналу можлива також оцінка помилок закруглення. За допомогою цього функціоналу тут отримані оцінки стійкості обчислень в маршевому методі розв'язування задачі про згин пластинки. Отримані результати можуть бути використані для оцінки часу обчислень в реальних ситуаціях в залежності від типу процесора і від варіанту програми, що використовується для обчислень.

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

В главі ІІ розглянуті сталі стійкості в метриках просторів L₂ та W₂¹ в класах трикутних сіток з рівномірно обмеженими кутами, а також з кутами, що прямують до нуля. В § 2.1 вказана роль сталих стійкості і побудований передобумовлювач. Він має просту структуру (блочно-трьохдіагональна матриця) і його можна економічно обертати. В § 2.2 і § 2.3 вивчені класи сіток з рівномірно обмеженими кутами в метриках просторів L₂, W₂¹ відповідно. В цих випадках наведені достатньо точні чисельні оцінки сталих. В § 2.4 вказані асимптотики поведінки границь відповідних квадратичних форм при вироджені трикутників.

Сформулюємо деякі з доведених в главі ІІ тверджень. Нехай - деяка область в R², розглянемо правильно тріангульований багатокутник . Перенумеруємо всі вершини нашої тріангуляції натуральними числами. Нехай _j = (x_j, y_j) - j-а вершина тріангуляції, j = 1, …, N. Множину трикутників, у яких j-а вершина загальна, назвемо барицентричною зіркою і позначимо Z_j. На Z_j задаємо кусково-лінійну і неперервну функцію _j(), що лінійна на кожному трикутнику, перетворюється в 1 в вершинах _j і в 0 в інших вершинах. В області будемо вважати цю функцію рівною 0. Розглянемо лінійну комбінацію курантовських функцій:

.

Очевидно:

де

Лема 1. Справедливі нерівності:

(8)

де mZ_j - площа багатокутника Z_j. Права оцінка в нерівності (8) є точною.

Для кожного трикутника тріангуляції введемо позначення: S_T - площа трикутнику, h₁ - висота, проведена з і-ої вершини на протилежну сторону, n₁ - одинична нормаль до сторони, протилежній вершині і, n_і = 1. z_і - звуження функції _і на трикутник Т,

(9)

При введених вище позначеннях інтеграли (9) обчислюються за формулами:

[z_і, z_і]= (S_T cos (n_і, n_j)) / (h_і, h_j) (10)

позначимо через ₀ сіткову функцію, що рівна 1 в усіх розглядуваних вузлах, ₀(_j) = 1, j=1, …, N.

Нехай V₀ - одновимірний підпростір сіткових функцій - сталих,

V₀ = {C₀: CR¹}.

Лема 2. Для сіткових функцій , що є ортогональні підпростору V₀, справедливі нерівності:

(11)

де (12)

(13)

Лема 3. Для довільних сіткових функцій справедливі нерівності:

(14)

де - ті ж сталі, що і в лемі 2, а сума означає підсумовування по всіx ребраx l_і,_j тріангуляції з кінцями у вершинах _і, _j.

Лема 4. Нехай визначаються формулою (13), , , - кути трикутника Т.

Якщо:

0, ₀, ₀, ₀ + ₀ = ₀ > 0 ₀ > 0,

то для власних чисел справедливі асимптотичні представлення:

= 3/4 + 0(³), ₂^Т = 1/ + 0(). (15)

Якщо = c, 0, c = const > 0, , то для вірні асимптотичні представлення:

= 3/4 с(с + 1)/(с ² + с+1) + 0(³), (16)

= (с ² + с + 1)/(с(с + 1)) 1/ + 0().

В главі ІІІ запропоновані різні способи побудови тріангуляції в залежності від гладкості границі. На підставі результатів другої глави тут отримані ефективні оцінки сталих стійкості для згаданих способів тріангуляції.

На площині R² вводиться тріангульована квадратна або прямокутна сітка і будується тріангуляція в приграничній смузі в залежності від гладкості границі. § 3.1 і § 3.2 дана ефективна оцінка знизу сталих при запропонованому способі тріангуляції в залежності від ширини вибраної приграничної смуги для двічі неперервно-диференційованої границі і від куту , що містить кусково-гладку границю. Чисельні результати наведені в таблицях 1-12 додатку 2. В § 3.3 пропонуються способи побудови тріангуляції для кусково-гладкої границі з внутрішніми нульовими кутами степеневих порядків. В § 3.4 запропонований алгоритм побудови - тріангуляції, в деякому сенсі близької до оптимальної.

Означення. Серед всіляких тріангуляцій, побудованих для заданої області, -тріангуляцією назвемо таку, для якої відношення max_1,2/ mіn_1,2 є мінімальним.

Сформулюємо деякі з доведених в главі ІІІ тверджень.

Як і раніше, нехай V₀ - одномірний підпростір сіткових функцій-сталих, - правильно тріангульований багатокутник.

Теорема 3. (Випадок двічі неперервно диференційованої границі). Нехай Г С ². Тоді для сіткових функцій , ортогональних до підпростору V₀, справедливі нерівності.

(17)

де ^(і), і = 1, …, 8 - сталі стійкості, що обчислені в дисертації (див. гл. ІІІ, § 3.1), а також в таблицях 1-4 додатку ІІ.

Теорема 4. (Кусково-гладка границя з ненульовим кутом).

Нехай Г належить до кусково-гладкої границі з ненульовим кутом, тоді для сіткових функцій , що ортогональні підпростору V₀, справедливі нерівності:

(18)

де ^(і)_1,2, і = 1, …., 8 - сталі стійкості, що наведені в главі ІІІ, § 3.2, а також в таблицях 5-12 додатку 2. Далі розглядається границя з нульовим кутом - зі ступеневими загостреннями вигляду y = х^a, a > 1.

Лема 5. Швидкість зменшення кутів нормальної тріангуляції в околі нульового кута ступеня а характеризується асимптоти кою:

= с h ^{1-1/ a} + …, (19)

де 2 а с ₀^{1-1/ a} с 2 а с ₁^{1-1/ a},

с ₀, с ₁ - деякі сталі.

Лема 6. Нехай Г - кусково-гладка границя з нульовим кутом, тоді справедливі нерівності:

₁(h) = 3/4 с h ^{1-1/ a} + 0(h ^{1-1/ a}),

₂(h) = 1/c h - ^{(1-1/ a)} - 3/4 с h ^{1-1/ a} + 0(h ^{1-1/ a}).

Теорема 5. (Кусково-гладка границя з нульовим кутом). Для сіткових функцій , що ортогональні підпростору V₀, справедливі нерівності:

(20)

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

Додаток ІІ містить таблиці сталих стійкості, отриманих в роботі, та малюнки. Обговорення чисельних результатів і структури таблиць міститься в § 3.1, § 3.2 глави ІІІ.

Додаток ІІІ містить ведення дійсності власних чисел квадратичних форм, що зустрічаються в главі ІІІ.

ОСНОВНІ РЕЗУЛЬТАТИ