Efektywna implementacja algorytmu macierzy tridiagonal

12

Rozwiązuję problem fizyczny za pomocą niejawnego schematu numerycznego. To prowadzi mnie do rozwiązania równania liniowego za pomocą macierzy tridiagonalnej. Algorytm kodowałem z Wikipedii. Zastanawiam się, czy istnieje wydajna biblioteka, która pozwala optymalnie rozwiązać tego rodzaju równanie. Ważną uwagą jest to, że sama matryca zmienia się tylko, gdy zmieniają się parametry systemu, więc miałem okazję wstępnie obliczyć kilka kroków algorytmu, aby uzyskać niezłą premię za wydajność. Używam C ++.

linear-algebra libraries c++ gmk
źródło
Jak duży jest system, czy musi być równoległy?
aterrel
1
Rozmiar zależy od wymaganej dokładności (od stu do dziesiątek tysięcy wartości). Teraz piszę na komputerze z jednym rdzeniem, ale można uzyskać dostęp do superkomputera uniwersyteckiego z wieloma dostępnymi procesorami, więc obsługa równoległości byłaby miła.
gmk

Odpowiedzi:

15

Prawdopodobnie powinieneś zacząć od implementacji LAPACK,? Gtsv, np . Dgtsv . Jeśli chcesz wersję z pamięcią rozproszoną, możesz zacząć od p? Gtsv ScaLAPACK.

EDYCJA: Ponieważ twoja macierz nie zmienia się zbyt często, możesz uniknąć zbędnego faktoryzowania macierzy tridiagonalnej, rozkładając procedurę LAPACK? Gtsv na etap faktoryzacji? Gttrf i etap rozwiązywania? Gttrs. W ScaLAPACK istnieją podobnie nazwane procedury, które służą temu samemu celowi.

Jack Poulson
źródło
Dziękuję, wygląda na to, czego potrzebuję. Spróbuję teraz uruchomić te procedury z mojego kodu.
gmk
1
Ponieważ wywołujesz go z C ++, zadeklaruj prototyp w zewnętrznym bloku „C” {}. W zależności od systemu może być konieczne dołączenie podkreślenia do nazwy rutyny.
Jack Poulson,
2

W przypadku rozproszonych systemów równoległych : Nie próbowałem ScaLAPACK, który ma równoległy trójdrożny solver, dla którego istnieją przykłady dostępne online . Z pewnym powodzeniem wypróbowałem metodę zaproponowaną przez Davida Moultona w publikacji LANL . Kodowanie tego może być więcej, niż chcesz, ale używając LAPACK, jest to proste.

Yann
źródło
1

Znalazłem interesujący algorytm rekurencyjny tutaj na stronie 975. Wygląda obiecująco, zastanawiam się, co mówią o tym bardziej doświadczeni ludzie.

Tiam
źródło
Przepisy numeryczne zawierają pewne błędy. Pod względem źródła kodów do użycia nie jest najlepszy, chociaż niektórzy uważają go za klasyczny. Byłbym zaskoczony, gdyby ScaLAPACK nie zaimplementował algorytmu co najmniej tak wydajnego jak rekurencyjna cykliczna redukcja.
Geoff Oxberry