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{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{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{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{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).

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.

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.