Dlaczego mój iteracyjny liniowy solver nie jest zbieżny?

26

Co może pójść nie tak, gdy zastosuje się wstępne metody Kryłowa z KSP ( pakiet solvera liniowego PETSc ) do rozwiązania rzadkiego układu liniowego, takiego jak te uzyskane przez dyskretyzację i linearyzację równań różniczkowych cząstkowych?

Jakie kroki mogę podjąć, aby ustalić, co jest nie tak z moim problemem?

Jakie zmiany mogę wprowadzić, aby skutecznie i skutecznie rozwiązać mój system liniowy?

Jed Brown
źródło
Czy to pytanie ma być pytaniem o iteracyjne solwery liniowe, konkretnie w PETSc (to, co bym zebrał z treści pytania), czy może być pytaniem o potencjalne awarie algorytmiczne iteracyjnych solverów liniowych w większości oprogramowania kontekst agnostyczny (co bym zebrał, patrząc na sam nagłówek)?
Geoff Oxberry
4
Ma petsctag. Metodologia jest ogólna, ale myślę, że odpowiedź byłaby mniej przydatna, gdyby każde „spróbuj tego” nie zawierało również „jak”. Alternatywnie „jak” musiałoby być znacznie dłuższe (i bardziej podatne na błędy dla widza), gdyby wymagało to wyjaśnienia w sposób niezależny od oprogramowania. Jeśli ktoś chce wyjaśnić, jak zrobić te wszystkie rzeczy przy użyciu innego pakietu, z przyjemnością zmienię pytanie w oprogramowanie i zmienię odpowiedź na stwierdzenie, że opisuje, co robić w PETSc. Uwaga: dodałem tę, która jest ulepszoną wersją FAQ, więc mogę polubić ludzi na tej stronie.
Jed Brown

Odpowiedzi:

26

Wstępna porada

  • Zawsze -ksp_converged_reason -ksp_monitor_true_residualstaraj się dowiedzieć, dlaczego metoda nie jest zbieżna.
  • Zmniejsz rozmiar problemu i liczbę procesów do jak najmniejszych, aby zademonstrować awarię. Często uzyskujesz wgląd, określając, jakie drobne problemy wykazują zachowanie, które powoduje awarię metody i skraca czas zawracania. Ponadto istnieją pewne techniki dochodzenia, które można zastosować tylko w przypadku małych systemów.
  • Jeśli problem pojawia się dopiero po wielu krokach czasowych, krokach kontynuacji lub krokach rozwiązywania nieliniowego, rozważ zapisanie stanu modelu, gdy wystąpi awaria, abyś mógł szybko eksperymentować.
  • Alternatywnie, szczególnie jeśli twoje oprogramowanie nie ma możliwości punktu kontrolnego, użyj -ksp_view_binarylub MatView()zapisz system liniowy, a następnie użyj kodu at, $PETSC_DIR/src/ksp/ksp/examples/tutorials/ex10.caby odczytać matrycę i rozwiązać ją (być może przy innej liczbie procesów). Wymaga to złożonej matrycy, więc jej użyteczność może być nieco ograniczona.
  • Istnieje wiele możliwych wyborów solvera (np. Nieskończona liczba dostępna w linii poleceń w PETSc ze względu na dowolną liczbę poziomów składu), zobacz to pytanie, aby uzyskać ogólne porady dotyczące wybierania solwerów liniowych.

Typowe przyczyny braku konwergencji KSP

  • Równania są przypadkowe osobliwe (np. Zapomniałem narzucić warunki brzegowe). Sprawdź to pod kątem małego problemu przy użyciu -pc_type svd -pc_svd_monitor. Wypróbuj również bezpośredni solver z -pc_type lu(za pośrednictwem pakietu innej firmy równolegle, np -pc_type lu -pc_factor_mat_solver_package superlu_dist.).
  • Równania są celowo pojedyncze (np. Stała pusta przestrzeń), ale metoda Kryłowa nie została poinformowana, patrz KSPSetNullSpace().
  • Równania są celowo pojedyncze i KSPSetNullSpace()zostały użyte, ale prawa strona nie jest spójna. Być może będziesz musiał zadzwonić MatNullSpaceRemove()po prawej stronie, zanim zadzwonisz KSPSolve().
  • Równania są nieokreślone, więc standardowe warunki wstępne nie działają. Zazwyczaj znasz to z fizyki, ale możesz to sprawdzić za pomocą -ksp_compute_eigenvalues -ksp_gmres_restart 1000 -pc_type none. W przypadku prostych problemów z siodłem spróbuj -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_detect_saddle_point. Aby uzyskać więcej informacji, zobacz Podręcznik użytkownika i stronę podręcznika PCFIELDSPLIT . W przypadku trudniejszych problemów zapoznaj się z literaturą, aby znaleźć solidne metody i zapytaj tutaj ( [email protected]lub [email protected]), czy chcesz uzyskać porady dotyczące ich wdrożenia. Na przykład zobacz to pytanie dla wysokiej częstotliwości Helmholtza. W przypadku skromnych rozmiarów problemów sprawdź, czy możesz żyć tylko za pomocą bezpośredniego solvera.
  • Jeśli metoda zbiega się w warstwie wstępnej kondycjonowanej, ale nie w reszcie resztkowej, to prawdopodobnie jednostka wstępna kondycjonująca jest pojedyncza lub prawie taka sama. Jest to powszechne w przypadku problemów z punktem siodłowym (np. Przepływ nieściśliwy) lub silnie niesymetrycznych operatorów (np. Problemy hiperboliczne o niskim Machu z dużymi krokami czasowymi).
  • Warunek wstępny jest zbyt słaby lub niestabilny. Sprawdź, czy -pc_type asm -sub_pc_type lupoprawia współczynnik konwergencji. Jeśli GMRES traci zbyt duży postęp w ponownym uruchomieniu, sprawdź, czy dłuższe ponowne uruchomienie nie pomaga -ksp_gmres_restart 300. Jeśli transpozycja jest dostępna, wypróbuj -ksp_type bcgsinne metody, które nie wymagają ponownego uruchomienia. (Należy zauważyć, że zbieżność z tymi metodami jest często błędna).
  • Matryca wstępnego kondycjonowania może nie znajdować się blisko (ewentualnie niezmontowanego) operatora. Spróbuj rozwiązać za pomocą solwera bezpośredniego, szeregowo -pc_type lulub równolegle, używając pakietu innej firmy (np. -pc_type lu -pc_factor_mat_solver_package superlu_distLub mumps). Metoda powinna zbiegać się w jednej iteracji, jeśli macierze są takie same, aw „małej” liczbie iteracji w przeciwnym razie. Spróbuj -snes_type testsprawdzić macierze, jeśli rozwiązujesz problem nieliniowy.
  • Warunek wstępny jest nieliniowy (np. Zagnieżdżone iteracyjne rozwiązanie), spróbuj -ksp_type fgmres or -ksp_type gcr.
  • Używasz geometrycznej wielosiatki, ale niektóre równania (często warunki brzegowe) nie są skalowane kompatybilnie między poziomami. Spróbuj -pc_mg_galerkinalgebraicznie skonstruować poprawnie wyskalowany operator zgrubny lub upewnij się, że wszystkie równania są skalowane w ten sam sposób, jeśli chcesz użyć ponownie zdyskretowanych poziomów zgrubnych.
  • Matryca jest bardzo źle uwarunkowana. Sprawdź numer warunku, stosując metody opisane tutaj . Spróbuj to poprawić, wybierając względne skalowanie komponentów / warunków brzegowych. Spróbować -ksp_diagonal_scale -ksp_diagonal_scale_fix. Być może zmień sformułowanie problemu, aby uzyskać bardziej przyjazne równania algebraiczne. Jeśli nie możesz poprawić skalowania, może być konieczne użycie bezpośredniego solvera.
  • Macierz jest nieliniowa (np. Oceniana przy użyciu różnicowania skończonego funkcji nieliniowej). Wypróbuj różne parametry różnicowania (np -mat_mffd_type ds.). Spróbuj użyć większą precyzję, aby różnicowanie dokładniejsze ./configure --with-precision=__float128 --download-f2cblaslapack. Sprawdź, czy zbiega się w „łatwiejszych” systemach parametrów.
  • W przypadku problemu niesymetrycznego stosowana jest metoda symetryczna.
  • Klasyczna Gram-Schmidt staje się niestabilna, spróbuj -ksp_gmres_modifiedgramschmidtlub zastosuj metodę, która inaczej ortogonalizuje, np -ksp_type gcr.
Jed Brown
źródło
16

Moja rada dla studentów polega na wypróbowaniu bezpośredniego rozwiązania w takich przypadkach. Powodem jest to, że istnieją dwie klasy powodów, dla których solver może się nie zbiegać: (i) matryca jest niepoprawna lub (ii) występuje problem z solverem / warunkiem wstępnym. Bezpośrednie solwery prawie zawsze dają coś, co można porównać z oczekiwanym rozwiązaniem, więc jeśli odpowiedź bezpośredniego solwera wygląda poprawnie, oznacza to, że problem dotyczy iteracyjnego solvera / warunku wstępnego. Z drugiej strony, jeśli odpowiedź wydaje się błędna, problem polega na złożeniu matrycy i prawej stronie.

Zwykle używam UMFPACK jako bezpośredniego solvera. Jestem pewien, że łatwo wypróbować coś podobnego z PETSC.

Wolfgang Bangerth
źródło
5
-pc_type lu -pc_factor_mat_solver_type umfpackużywać UMFPACK (lub w -pc_type cholesky -pc_factor_mat_solver_package cholmodprzypadku problemów z SPD) za pośrednictwem PETSc, ale pamiętaj, że UMFPACK i CHOLMOD są szeregowe. Dla równoległej stosowania -pc_factor_mat_solver_package superlu_distlub mumps, pastix, spooles.
Jed Brown
2
Dla jasności kompletny zestaw opcji do używania (np.) superlu_distByłby -ksp_type preonly -pc_type lu -pc_factor_mat_solver_package superlu_dist. Czy to prawda?
Leon Avery
Nie wiem Co się stanie, jeśli to zrobisz?
Wolfgang Bangerth,