Jak skutecznie zaimplementować warunki brzegowe Dirichleta w globalnych macierzach sztywnych elementów skończonych

Zastanawiam się, w jaki sposób warunki brzegowe Dirichleta w globalnych macierzach rzadkich elementów skończonych są faktycznie skutecznie wdrażane. Powiedzmy na przykład, że nasza globalna macierz elementów skończonych była:

K = [\begin{matrix} 5 & 2 & 0 & - 1 & 0 \\ 2 & 4 & 1 & 0 & 0 \\ 0 & 1 & 6 & 3 & 2 \\ - 1 & 0 & 3 & 7 & 0 \\ 0 & 0 & 2 & 0 & 3 \end{matrix}] and right-hand side vector b = [\begin{matrix} b 1 \\ b 2 \\ b 3 \\ b 4 \\ b 5 \end{matrix}]

$K = \begin{bmatrix} 5 & 2 & 0 & -1 & 0 \\ 2 & 4 & 1 & 0 & 0 \\ 0 & 1 & 6 & 3 & 2 \\ -1 & 0 & 3 & 7 & 0 \\ 0 & 0 & 2 & 0 & 3 \end{bmatrix} \hspace{5mm}\text{and right-hand side vector}\hspace{5mm} b = \begin{bmatrix} b1 \\ b2 \\ b3 \\ b4 \\ b5 \\ \end{bmatrix}$

Następnie, aby zastosować warunek Dirichleta w pierwszym węźle ( $x_{1}=c$ ) wyzerowalibyśmy pierwszy rząd, postawiliśmy 1 na $K_{11}$ i odejmij pierwszą kolumnę od prawej strony. Na przykład nasz system stałby się:

K = [\begin{matrix} 1 & 0 & 0 & 0 & 0 \\ 0 & 4 & 1 & 0 & 0 \\ 0 & 1 & 6 & 3 & 2 \\ 0 & 0 & 3 & 7 & 0 \\ 0 & 0 & 2 & 0 & 3 \end{matrix}] and right-hand side vector b = [\begin{matrix} c \\ b 2 - 2 \times c \\ b 3 - 0 \times c \\ b 4 + 1 \times c \\ b 5 - 0 \times c \end{matrix}]

$K = \begin{bmatrix} 1 & 0 & 0 & 0 & 0 \\ 0 & 4 & 1 & 0 & 0 \\ 0 & 1 & 6 & 3 & 2 \\ 0 & 0 & 3 & 7 & 0 \\ 0 & 0 & 2 & 0 & 3 \end{bmatrix} \hspace{5mm}\text{and right-hand side vector}\hspace{5mm} b = \begin{bmatrix} c \\ b2-2\times{c} \\ b3-0\times{c} \\ b4+1\times{c} \\ b5-0\times{c} \\ \end{bmatrix}$

W teorii wszystko jest dobrze i dobrze, ale jeśli nasza macierz K jest przechowywana w formacie skompresowanego wiersza (CRS), wówczas przenoszenie kolumn po prawej stronie staje się kosztowne w przypadku dużych systemów (wiele węzłów jest dirichlet). Alternatywą byłoby nie przesunięcie kolumn odpowiadających warunkowi Dirichleta na prawą stronę, tzn. Nasz system stałby się:

K = [\begin{matrix} 1 & 0 & 0 & 0 & 0 \\ 2 & 4 & 1 & 0 & 0 \\ 0 & 1 & 6 & 3 & 2 \\ - 1 & 0 & 3 & 7 & 0 \\ 0 & 0 & 2 & 0 & 3 \end{matrix}] and right-hand side vector b = [\begin{matrix} c \\ b 2 \\ b 3 \\ b 4 \\ b 5 \end{matrix}]

$K = \begin{bmatrix} 1 & 0 & 0 & 0 & 0 \\ 2 & 4 & 1 & 0 & 0 \\ 0 & 1 & 6 & 3 & 2 \\ -1 & 0 & 3 & 7 & 0 \\ 0 & 0 & 2 & 0 & 3 \end{bmatrix} \hspace{5mm}\text{and right-hand side vector}\hspace{5mm} b = \begin{bmatrix} c \\ b2 \\ b3 \\ b4 \\ b5 \\ \end{bmatrix}$

Ma to jednak poważną wadę polegającą na tym, że system nie jest już symetryczny i dlatego nie możemy już używać wstępnie przygotowanego gradientu sprzężonego (lub innych symetrycznych solverów). Jednym z interesujących rozwiązań, na jakie natrafiłem, jest „Metoda dużych liczb”, którą znalazłem w książce „Programowanie elementów skończonych w Javie” Gennadiya Nikishkova. Ta metoda wykorzystuje fakt, że podwójna precyzja zawiera jedynie około 16 cyfr dokładności. Zamiast wstawiać 1 w $K_{11}$ pozycja umieszczamy dużą liczbę. Na przykład nasz system staje się:

K = [\begin{matrix} 1.0 e 64 & 2 & 0 & - 1 & 0 \\ 2 & 4 & 1 & 0 & 0 \\ 0 & 1 & 6 & 3 & 2 \\ - 1 & 0 & 3 & 7 & 0 \\ 0 & 0 & 2 & 0 & 3 \end{matrix}] and right-hand side vector b = [\begin{matrix} c \times 1.0 e 64 \\ b 2 \\ b 3 \\ b 4 \\ b 5 \end{matrix}]

$K = \begin{bmatrix} 1.0e64 & 2 & 0 & -1 & 0 \\ 2 & 4 & 1 & 0 & 0 \\ 0 & 1 & 6 & 3 & 2 \\ -1 & 0 & 3 & 7 & 0 \\ 0 & 0 & 2 & 0 & 3 \end{bmatrix} \hspace{5mm}\text{and right-hand side vector}\hspace{5mm} b = \begin{bmatrix} c\times{1.0e64} \\ b2 \\ b3 \\ b4 \\ b5 \\ \end{bmatrix}$

Zaletą tej metody jest to, że utrzymuje symetrię matrycy, a jednocześnie jest bardzo wydajna w przypadku rzadkich formatów pamięci. Moje pytania są zatem następujące:

W jaki sposób warunki brzegowe Dirichleta są zazwyczaj realizowane w kodach elementów skończonych dla ciepła / płynów? Czy ludzie zwykle stosują metodę dużych liczb, czy robią coś innego? Czy jest jakaś wada metody dużych liczb, którą ktoś może zobaczyć? Zakładam, że w większości komercyjnych i niekomercyjnych kodów istnieje prawdopodobnie jakaś standardowa wydajna metoda, która rozwiązuje ten problem (oczywiście nie oczekuję, że ludzie będą znali wszystkie wewnętrzne działania każdego komercyjnego rozwiązania elementu skończonego, ale problem ten wydaje się podstawowy / fundamentalny wystarczy, że ktoś prawdopodobnie pracował nad takimi projektami i mógł udzielić wskazówek).

finite-element sparse-matrix boundary-conditions

— James
źródło

Czy masz dowody na to, że znacznie Cię to spowalnia?

— Bill Barth

@BillBarth Tak, chociaż zawsze istnieje szansa, że robię coś nieefektywnie. Sam Gennadily pisze, że chociaż metoda jawna jest łatwa dla pełnych tablic 2D, „.. nie zawsze jest łatwy dostęp do wierszy i kolumn macierzy, gdy macierz ma kompaktowy format”. sugerując, że precyzyjna metoda może być trudniejsza do wydajnego wdrożenia. Ponieważ mój kod jest obecnie napisany, jawna metoda może zająć więcej czasu niż samo rozwiązanie.

— James

zrób to tak, jak mówi Wolfgang i zastosuj warunki brzegowe do matryc elementów przed złożeniem.

— Bill Barth,

@BillBarth Tak, myślę, że to zrobię. Jego filmy są niesamowite! Właśnie zostawiłem mu komentarz / pytanie dotyczące tego, czy musisz składać globalne macierze za każdym razem, po czym myślę, że zaakceptuję jego odpowiedź.

— James

W deal.II ( http://www.dealii.org - wyłączenie odpowiedzialności: Jestem jednym z głównych autorów tej biblioteki), eliminujemy całe wiersze i kolumny i ogólnie nie jest to zbyt drogie. Sztuką jest wykorzystanie faktu, że wzorzec rzadkości jest zwykle symetryczny, więc wiesz, w których rzędach musisz sprawdzić, eliminując całą kolumnę.

Moim zdaniem lepsze podejście polega na wyeliminowaniu tych wierszy i kolumn w macierzach komórek, zanim zostaną one dodane do macierzy globalnej. Tam pracujesz z pełnymi matrycami, więc wszystko jest wydajne.

Nigdy nie słyszałem o podejściu do wielkich liczb i nie zastosowałbym go, ponieważ z pewnością doprowadzi to do bardzo źle uwarunkowanych problemów.

Dla porównania, algorytmy, których używamy w deal.II, są opisane koncepcyjnie w wykładach 21.6 i 21.65 na stronie http://www.math.colostate.edu/~bangerth/videos.html . Dokładnie pasują do twojego opisu.

— Wolfgang Bangerth
źródło

Czy w przypadku problemu zależnego od czasu (powiedzmy równania ciepła) ponownie składasz macierz globalną za każdym razem? Powód, dla którego pytam, jest taki, że w przypadku niezerowych warunków Dirichleta potrzebujesz informacji z oryginalnej globalnej macierzy podczas modyfikowania prawej strony, ale jeśli wyzerowałeś te kolumny podczas poprzedniego pomiaru czasu, informacje te zostaną utracone (chyba że je zapiszesz w dodatkowych tablicach). Nie stanowiłoby to problemu, gdyby globalna macierz była ponownie składana przy każdym kroku czasowym, co właśnie rozważam i co i tak musiałoby zostać zrobione, gdyby zastosowano siatkę adaptacyjną.

— James

To zależy od zastosowania. Wszystkie „duże” kody rozwiązują nieliniowe problemy zależne od czasu, a dla nich jasne jest, że musisz ponownie złożyć jeden lub drugi sposób. W przypadku kodów liniowych możesz po prostu zapisać oryginalną matrycę i za każdym razem kopiować ją gdzie indziej, zastosować warunki brzegowe, a następnie użyć jej w solwerze. To po prostu wymaga więcej pamięci, ale poza tym jest tanie.

— Wolfgang Bangerth,

Ach, widzę, że tak podejrzewałem. Wdrożę, jak zasugerowałeś. Ok, to za twoją pomoc. Te filmy instruktażowe deallii są naprawdę dobre!

— James

Zero BC są łatwe. W przypadku niezerowych BC można także użyć mnożników Lagrange'a. Np. Patrz tutaj . Jedną z zalet LM jest to, że możesz użyć dowolnego równania ograniczenia, chociaż system staje się nieokreślony, więc potrzebujesz odpowiedniego rozwiązania.

— stali
źródło