W jaki sposób BLAS osiąga tak ekstremalną wydajność?

108

Z ciekawości zdecydowałem się porównać moją własną funkcję mnożenia macierzy z implementacją BLAS-a ... Wynik był, powiem, najmniej zaskoczony:

Implementacja niestandardowa, 10 prób mnożenia macierzy 1000x1000:

Took: 15.76542 seconds.

Implementacja BLAS, 10 prób mnożenia macierzy 1000x1000:

Took: 1.32432 seconds.

To jest używanie liczb zmiennoprzecinkowych o pojedynczej precyzji.

Moje wdrożenie:

template<class ValT>
void mmult(const ValT* A, int ADim1, int ADim2, const ValT* B, int BDim1, int BDim2, ValT* C)
{
    if ( ADim2!=BDim1 )
        throw std::runtime_error("Error sizes off");

    memset((void*)C,0,sizeof(ValT)*ADim1*BDim2);
    int cc2,cc1,cr1;
    for ( cc2=0 ; cc2<BDim2 ; ++cc2 )
        for ( cc1=0 ; cc1<ADim2 ; ++cc1 )
            for ( cr1=0 ; cr1<ADim1 ; ++cr1 )
                C[cc2*ADim2+cr1] += A[cc1*ADim1+cr1]*B[cc2*BDim1+cc1];
}

Mam dwa pytania:

  1. Biorąc pod uwagę, że mnożenie macierzy-macierzy mówi: nxm * mxn wymaga n * n * m mnożenia, czyli w przypadku powyżej 1000 ^ 3 lub 1e9 operacji. W jaki sposób na moim procesorze 2,6 GHz BLAS może wykonać 10 * 1e9 operacji w 1,32 sekundy? Nawet jeśli mnożenie było pojedynczą operacją i nic innego nie zostało zrobione, powinno to zająć ~ 4 sekundy.
  2. Dlaczego moja implementacja jest o wiele wolniejsza?
DeusAduro
źródło
17
BLAS został zoptymalizowany z jednej strony iz drugiej przez specjalistów w tej dziedzinie. Zakładam, że korzysta z jednostki zmiennoprzecinkowej SIMD na twoim chipie i robi wiele sztuczek, aby poprawić zachowanie buforowania ...
dmckee --- kociak ex-moderator
3
Jak jednak wykonać operacje 1E10 na procesorze 2,63E9 na sekundę w 1,3 sekundy?
DeusAduro
9
Wiele jednostek wykonawczych, układanie rur i pojedyncze rozkazowe dane wielokrotne ((SIMD), co oznacza wykonywanie tej samej operacji na więcej niż jednej parze operandów w tym samym czasie). Niektóre kompilatory mogą kierować jednostki SIMD na zwykłe chipy, ale prawie zawsze musisz jawnie włączyć, a warto wiedzieć, jak to wszystko działa ( en.wikipedia.org/wiki/SIMD ). Ubezpieczenie od utraty pamięci podręcznej jest prawie na pewno najtrudniejszą częścią.
dmckee --- kociak ex-moderator
13
Założenie jest błędne. Znane są lepsze algorytmy, patrz Wikipedia.
MSalters
2
@DeusAduro: W mojej odpowiedzi na pytanie Jak napisać produkt macierzowy, który może konkurować z firmą Eigen? Opublikowałem mały przykład, jak zaimplementować produkt macierzowo-macierzowy o wydajności pamięci podręcznej.
Michael Lehn

Odpowiedzi:

141

Dobrym punktem wyjścia jest wspaniała książka „ The Science of Programming Matrix Computations” autorstwa Roberta A. van de Geijna i Enrique S. Quintana-Ortí. Zapewniają bezpłatną wersję do pobrania.

BLAS jest podzielony na trzy poziomy:

  • Poziom 1 definiuje zestaw funkcji algebry liniowej, które działają tylko na wektorach. Funkcje te czerpią korzyści z wektoryzacji (np. Przy użyciu SSE).

  • Funkcje poziomu 2 to operacje macierzowo-wektorowe, np. Jakiś iloczyn macierzowo-wektorowy. Funkcje te można by zaimplementować w kategoriach funkcji poziomu 1. Możesz jednak zwiększyć wydajność tej funkcji, jeśli możesz zapewnić dedykowaną implementację, która wykorzystuje architekturę wieloprocesorową ze współdzieloną pamięcią.

  • Funkcje poziomu 3 to operacje takie jak iloczyn macierzy-macierzy. Ponownie możesz zaimplementować je w kategoriach funkcji Level2. Ale funkcje Level3 wykonują operacje O (N ^ 3) na danych O (N ^ 2). Jeśli więc Twoja platforma ma hierarchię pamięci podręcznej, możesz zwiększyć wydajność, jeśli zapewnisz dedykowaną implementację, która jest zoptymalizowana pod kątem pamięci podręcznej / przyjazna dla pamięci podręcznej . Jest to ładnie opisane w książce. Głównym wzmocnieniem funkcji Level3 jest optymalizacja pamięci podręcznej. To przyspieszenie znacznie przewyższa drugie wzmocnienie z równoległości i innych optymalizacji sprzętu.

Nawiasem mówiąc, większość (lub nawet wszystkie) wysokowydajnych implementacji BLAS NIE jest zaimplementowanych w Fortranie. ATLAS jest zaimplementowany w C. GotoBLAS / OpenBLAS jest zaimplementowany w C, a jego części krytyczne dla wydajności w Assemblerze. W Fortranie zaimplementowano tylko referencyjną implementację BLAS. Jednak wszystkie te implementacje BLAS zapewniają interfejs Fortran w taki sposób, że można go połączyć z LAPACK (LAPACK zyskuje całą swoją wydajność z BLAS).

Zoptymalizowane kompilatory odgrywają pod tym względem niewielką rolę (a dla GotoBLAS / OpenBLAS kompilator nie ma żadnego znaczenia).

IMHO żadna implementacja BLAS nie wykorzystuje algorytmów, takich jak algorytm Coppersmith-Winograd lub algorytm Strassena. Nie jestem do końca pewien, dlaczego tak jest, ale zgaduję:

  • Może nie jest możliwe zapewnienie implementacji zoptymalizowanej pod kątem pamięci podręcznej tych algorytmów (tj. Straciłbyś więcej niż wygrał)
  • Te algorytmy nie są stabilne numerycznie. Ponieważ BLAS jest obliczeniowym jądrem LAPACK, nie można tego zrobić.

Edycja / aktualizacja:

Nowym i przełomowym dokumentem na ten temat są dokumenty BLIS . Są wyjątkowo dobrze napisane. Na wykładzie "Podstawy oprogramowania do obliczeń o wysokiej wydajności" zaimplementowałem produkt macierzowo-macierzowy po ich artykule. Właściwie zaimplementowałem kilka wariantów produktu macierzowo-macierzowego. Najprostsze warianty są w całości napisane w czystym C i mają mniej niż 450 linii kodu. Wszystkie inne warianty tylko optymalizują pętle

    for (l=0; l<MR*NR; ++l) {
        AB[l] = 0;
    }
    for (l=0; l<kc; ++l) {
        for (j=0; j<NR; ++j) {
            for (i=0; i<MR; ++i) {
                AB[i+j*MR] += A[i]*B[j];
            }
        }
        A += MR;
        B += NR;
    }

Ogólna wydajność produktu macierz-macierz zależy tylko od tych pętli. Spędza się tu około 99,9% czasu. W pozostałych wariantach użyłem elementów wewnętrznych i kodu asemblera, aby poprawić wydajność. Możesz zobaczyć samouczek przedstawiający wszystkie warianty tutaj:

ulmBLAS: Samouczek dotyczący GEMM (produkt Matrix-Matrix)

Wraz z dokumentami BLIS dość łatwo można zrozumieć, w jaki sposób biblioteki takie jak Intel MKL mogą uzyskać taką wydajność. I dlaczego nie ma znaczenia, czy używasz pamięci głównej w postaci wierszy czy kolumn!

Końcowe testy porównawcze są tutaj (nazwaliśmy nasz projekt ulmBLAS):

Benchmarki dla ulmBLAS, BLIS, MKL, openBLAS i Eigen

Kolejna edycja / aktualizacja:

Napisałem również kilka poradników na temat tego, jak BLAS jest używany do rozwiązywania problemów z algebry liniowej numerycznej, takich jak rozwiązywanie układu równań liniowych:

Wysokowydajna faktoryzacja LU

(Ta faktoryzacja LU jest na przykład używana przez Matlab do rozwiązywania układu równań liniowych).

Mam nadzieję, że znajdę czas na rozszerzenie tego samouczka, aby opisać i zademonstrować, jak zrealizować wysoce skalowalną równoległą implementację faktoryzacji LU, jak w PLASMA .

OK, gotowe: kodowanie równoległej faktoryzacji LU zoptymalizowanej pod kątem pamięci podręcznej

PS: Zrobiłem też kilka eksperymentów nad poprawą wydajności uBLAS. W rzeczywistości jest całkiem proste, aby zwiększyć (tak, grać słowami :)) wydajność uBLAS:

Eksperymenty na uBLAS .

Tutaj podobny projekt z BLAZE :

Eksperymenty na BLAZE .

Michael Lehn
źródło
3
Nowy link do „Benchmarki dla ulmBLAS, BLIS, MKL, openBLAS i Eigen”: apfel.mathematik.uni-ulm.de/~lehn/ulmBLAS/#toc3
Ahmed Fasih
Okazuje się, że ESSL firmy IBM wykorzystuje odmianę algorytmu Strassena - ibm.com/support/knowledgecenter/en/SSFHY8/essl_welcome.html
ben-albrecht
2
większość linków jest martwa
Aurélien Pierre
PDF z TSoPMC można znaleźć na stronie autora, pod adresem cs.utexas.edu/users/rvdg/tmp/TSoPMC.pdf
Alex Shpilkin
Chociaż algorytm Coppersmith-Winograd ma ładną złożoność czasową na papierze, notacja Big O ukrywa bardzo dużą stałą, więc zaczyna być opłacalna tylko dla śmiesznie dużych macierzy.
DiehardTheTryhard
26

Przede wszystkim BLAS to tylko interfejs z około 50 funkcjami. Istnieje wiele konkurencyjnych implementacji interfejsu.

Najpierw wspomnę o rzeczach, które są w dużej mierze niezwiązane:

  • Fortran vs C, nie ma znaczenia
  • Zaawansowane algorytmy macierzowe, takie jak Strassen, nie używają ich, ponieważ nie pomagają w praktyce

Większość implementacji dzieli każdą operację na małe macierze lub operacje wektorowe w mniej lub bardziej oczywisty sposób. Na przykład duże mnożenie macierzy 1000x1000 może zostać podzielone na sekwencję mnożenia macierzy 50x50.

Te operacje o małych wymiarach o stałym rozmiarze (zwane jądrem) są zakodowane na stałe w kodzie asemblera specyficznym dla procesora przy użyciu kilku funkcji procesora docelowego:

  • Instrukcje w stylu SIMD
  • Równoległość poziomu instrukcji
  • Świadomość pamięci podręcznej

Co więcej, jądra te mogą być wykonywane równolegle względem siebie przy użyciu wielu wątków (rdzeni procesora), w typowym wzorcu projektowym zmniejszania map.

Spójrz na ATLAS, który jest najczęściej używaną implementacją BLASa typu open source. Ma wiele różnych konkurujących ze sobą jąder, a podczas procesu budowania biblioteki ATLAS uruchamia konkurencję między nimi (niektóre są nawet sparametryzowane, więc to samo jądro może mieć różne ustawienia). Próbuje różnych konfiguracji, a następnie wybiera najlepszą dla konkretnego systemu docelowego.

(Wskazówka: dlatego, jeśli używasz ATLAS, lepiej jest zbudować i dostroić bibliotekę ręcznie dla konkretnej maszyny, niż używając wstępnie zbudowanej).

Andrew Tomazos
źródło
ATLAS nie jest już najczęściej używaną implementacją BLAS typu open source. Został przebity przez OpenBLAS (rozwidlenie GotoBLAS) i BLIS (refaktoryzacja GotoBLAS).
Robert van de Geijn
1
@ ulaff.net: To może. To zostało napisane 6 lat temu. Myślę, że obecnie najszybszą implementacją BLAS (oczywiście na Intelu) jest Intel MKL, ale nie jest to oprogramowanie typu open source.
Andrew Tomazos,
14

Po pierwsze, istnieją bardziej wydajne algorytmy mnożenia macierzy niż ten, którego używasz.

Po drugie, twój procesor może wykonać więcej niż jedną instrukcję naraz.

Twój procesor wykonuje 3-4 instrukcje na cykl, a jeśli używane są jednostki SIMD, każda instrukcja przetwarza 4 liczby zmiennoprzecinkowe lub 2 podwójne. (oczywiście ta liczba również nie jest dokładna, ponieważ procesor zwykle może przetwarzać tylko jedną instrukcję SIMD na cykl)

Po trzecie, twój kod jest daleki od optymalnego:

  • Używasz surowych wskaźników, co oznacza, że ​​kompilator musi założyć, że mogą one alias. Istnieją słowa kluczowe lub flagi specyficzne dla kompilatora, które można określić, aby poinformować kompilator, że nie mają aliasów. Alternatywnie powinieneś użyć innych typów niż surowe wskaźniki, które rozwiązują problem.
  • Niszczysz pamięć podręczną, wykonując naiwne przeglądanie każdego wiersza / kolumny macierzy wejściowych. Możesz użyć blokowania, aby wykonać jak najwięcej pracy na mniejszym bloku macierzy, który mieści się w pamięci podręcznej procesora, przed przejściem do następnego bloku.
  • W przypadku zadań czysto numerycznych Fortran jest prawie nie do pokonania, a C ++ wymaga wiele wysiłku, aby uzyskać podobną prędkość. Można to zrobić i jest kilka bibliotek demonstrujących to (zwykle przy użyciu szablonów wyrażeń), ale nie jest to trywialne i nie dzieje się tak po prostu .
jalf
źródło
Dzięki, dodałem poprawny kod ograniczania zgodnie z sugestią Justicle, nie widziałem dużej poprawy, podoba mi się pomysł blokowy. Z ciekawości, bez znajomości rozmiaru pamięci podręcznej procesora, jak można dobrać optymalny kod?
DeusAduro
2
Ty nie. Aby uzyskać optymalny kod, musisz znać rozmiar pamięci podręcznej procesora. Oczywiście wadą jest to, że skutecznie kodujesz swój kod na sztywno, aby uzyskać najlepszą wydajność na jednej rodzinie procesorów.
jalf
2
Przynajmniej pętla wewnętrzna pozwala uniknąć obciążeń krokowych. Wygląda na to, że jest to napisane dla jednej macierzy, która jest już transponowana. Dlatego jest „tylko” o jeden rząd wielkości wolniejszy niż BLAS! Ale tak, wciąż się rzuca z powodu braku blokowania pamięci podręcznej. Czy na pewno Fortran bardzo by pomógł? Myślę, że wszystko, co można tutaj zyskać, to to, że restrict(brak aliasingu) jest wartością domyślną, w przeciwieństwie do C / C ++. (I niestety ISO C ++ nie ma restrictsłowa kluczowego, więc musisz go używać __restrict__na kompilatorach, które zapewniają to jako rozszerzenie).
Peter Cordes
11

Nie wiem konkretnie o implementacji BLAS-a, ale istnieją bardziej wydajne alogorytmy dla mnożenia macierzy, które mają lepszą złożoność niż O (n3). Dobrze znanym jest Algorytm Strassena

softveda
źródło
8
Algorytm Strassena nie jest używany w numeryce z dwóch powodów: 1) nie jest stabilny. 2) Oszczędzasz niektóre obliczenia, ale wiąże się to z ceną, jaką możesz wykorzystać hierarchie pamięci podręcznej. W praktyce tracisz nawet wydajność.
Michael Lehn
4
W celu praktycznej implementacji algorytmu Strassen, ściśle opartego na kodzie źródłowym biblioteki BLAS, opublikowano niedawno publikację: „ Strassen Algorithm Reloaded ” w SC16, która osiąga wyższą wydajność niż BLAS, nawet dla problemu o rozmiarze 1000x1000.
Jianyu Huang
4

Większość argumentów na drugie pytanie - asembler, dzielenie na bloki itp. (Ale nie mniej niż algorytmy N ^ 3, są one naprawdę nadmiernie rozbudowane) - odgrywa rolę. Ale niska prędkość twojego algorytmu jest spowodowana głównie rozmiarem macierzy i niefortunnym ułożeniem trzech zagnieżdżonych pętli. Twoje macierze są tak duże, że nie mieszczą się od razu w pamięci podręcznej. Możesz przestawić pętle w taki sposób, aby jak najwięcej było zrobionych w wierszu w pamięci podręcznej, w ten sposób radykalnie zmniejszając odświeżanie pamięci podręcznej (przy okazji dzielenie na małe bloki ma efekt analogowy, najlepiej, jeśli pętle nad blokami są ułożone podobnie). Następuje modelowa implementacja macierzy kwadratowych. Na moim komputerze jego zużycie czasu wyniosło około 1:10 w porównaniu ze standardową implementacją (taką jak Twoja). Innymi słowy: nigdy nie programuj mnożenia macierzy wzdłuż "

    void vector(int m, double ** a, double ** b, double ** c) {
      int i, j, k;
      for (i=0; i<m; i++) {
        double * ci = c[i];
        for (k=0; k<m; k++) ci[k] = 0.;
        for (j=0; j<m; j++) {
          double aij = a[i][j];
          double * bj = b[j];
          for (k=0; k<m; k++)  ci[k] += aij*bj[k];
        }
      }
    }

Jeszcze jedna uwaga: ta implementacja jest nawet lepsza na moim komputerze niż zastąpienie wszystkich przez procedurę BLAS cblas_dgemm (wypróbuj ją na swoim komputerze!). Ale znacznie szybsze (1: 4) jest bezpośrednie wywołanie dgemm_ z biblioteki Fortran. Myślę, że ta procedura nie jest w rzeczywistości Fortranem, ale kodem assemblera (nie wiem, co jest w bibliotece, nie mam źródeł). Zupełnie niejasne jest dla mnie, dlaczego cblas_dgemm nie jest tak szybki, skoro według mojej wiedzy jest tylko opakowaniem dla dgemm_.

Wolfgang Jansen
źródło
3

To realistyczne przyspieszenie. Aby zapoznać się z przykładem tego, co można zrobić za pomocą asemblera SIMD w kodzie C ++, zobacz kilka przykładowych funkcji macierzy iPhone'a - były one ponad 8 razy szybsze niż wersja C i nie są nawet „zoptymalizowane” montażem - nie ma jeszcze wykładania rur i tam to niepotrzebne operacje na stosie.

Również twój kod nie jest „ ogranicz poprawny ” - skąd kompilator wie, że modyfikując C, nie modyfikuje A i B?

Justicle
źródło
Jasne, jeśli nazwałaś funkcję jak mmult (A ..., A ..., A); na pewno nie uzyskasz oczekiwanego rezultatu. Znowu jednak nie próbowałem pokonać / ponownie zaimplementować BLAS-a, po prostu sprawdzając, jak naprawdę jest szybki, więc nie chodziło o sprawdzanie błędów, tylko o podstawową funkcjonalność.
DeusAduro
3
Przepraszam, żeby było jasne, mówię, że jeśli ustawisz „ogranicz” na swoich wskaźnikach, otrzymasz znacznie szybszy kod. Dzieje się tak, ponieważ za każdym razem, gdy modyfikujesz C, kompilator nie musi przeładowywać A i B - radykalnie przyspieszając wewnętrzną pętlę. Jeśli mi nie wierzysz, sprawdź demontaż.
Justicle
@DeusAduro: To nie jest sprawdzanie błędów - możliwe, że kompilator nie może zoptymalizować dostępu do tablicy B [] w pętli wewnętrznej, ponieważ może nie być w stanie dowiedzieć się, że wskaźniki A i C nigdy nie są aliasami B szyk. Gdyby istniał aliasowanie, wartość w tablicy B mogłaby ulec zmianie podczas wykonywania pętli wewnętrznej. Wyciągnięcie dostępu do wartości B [] z pętli wewnętrznej i umieszczenie jej w zmiennej lokalnej może umożliwić kompilatorowi uniknięcie ciągłego dostępu do B [].
Michael Burr
1
Hmmm, więc najpierw spróbowałem użyć słowa kluczowego „__restrict” w VS 2008, zastosowanego do A, B i C. Wynik nie wykazał żadnej zmiany. Jednak przeniesienie dostępu do B, z najbardziej wewnętrznej pętli do pętli na zewnątrz, poprawiło czas o ~ 10%.
DeusAduro
1
Przepraszam, nie jestem pewien co do VC, ale w przypadku GCC musisz włączyć -fstrict-aliasing. Jest też lepsze wyjaśnienie „ograniczenia” tutaj: cellperformance.beyond3d.com/articles/2006/05/…
Justicle
2

W odniesieniu do oryginalnego kodu w mnożeniu MM, odwołanie do pamięci dla większości operacji jest główną przyczyną złej wydajności. Pamięć działa 100-1000 razy wolniej niż pamięć podręczna.

Większość przyspieszenia wynika z zastosowania technik optymalizacji pętli dla tej funkcji potrójnej pętli w mnożeniu MM. Stosowane są dwie główne techniki optymalizacji pętli; rozwijanie i blokowanie. Jeśli chodzi o rozwijanie, rozwijamy dwie zewnętrzne najbardziej pętle i blokujemy je w celu ponownego wykorzystania danych w pamięci podręcznej. Odwijanie pętli zewnętrznej pomaga tymczasowo zoptymalizować dostęp do danych, zmniejszając liczbę odwołań do pamięci do tych samych danych w różnym czasie podczas całej operacji. Zablokowanie indeksu pętli pod określonym numerem pomaga w zachowaniu danych w pamięci podręcznej. Możesz wybrać optymalizację dla pamięci podręcznej L2 lub L3.

https://en.wikipedia.org/wiki/Loop_nest_optimization

Pari Rajaram
źródło
-24

Z wielu powodów.

Po pierwsze, kompilatory Fortrana są wysoce zoptymalizowane, a język pozwala im takimi być. C i C ++ są bardzo luźne pod względem obsługi tablic (np. Przypadek wskaźników odnoszących się do tego samego obszaru pamięci). Oznacza to, że kompilator nie może z góry wiedzieć, co robić i jest zmuszony do utworzenia kodu ogólnego. W Fortranie sprawy są bardziej uproszczone, a kompilator ma lepszą kontrolę nad tym, co się dzieje, co pozwala mu na większą optymalizację (np. Przy użyciu rejestrów).

Inną rzeczą jest to, że Fortran przechowuje dane kolumnowo, podczas gdy C przechowuje dane wierszowo. Nie sprawdziłem kodu, ale uważaj na sposób wykonywania produktu. W C musisz skanować mądrze wierszami: w ten sposób skanujesz swoją tablicę wzdłuż ciągłej pamięci, zmniejszając błędy pamięci podręcznej. Brak pamięci podręcznej jest pierwszym źródłem nieefektywności.

Po trzecie, zależy to od używanej implementacji blas. Niektóre implementacje mogą być napisane w asemblerze i zoptymalizowane dla konkretnego używanego procesora. Wersja netlib jest napisana w Fortran 77.

Ponadto wykonujesz wiele operacji, większość z nich jest powtarzana i zbędna. Wszystkie te mnożenia w celu uzyskania indeksu są szkodliwe dla wydajności. Naprawdę nie wiem, jak to się robi w BLAS, ale jest wiele sztuczek, aby zapobiec kosztownym operacjom.

Na przykład możesz w ten sposób przerobić swój kod

template<class ValT>
void mmult(const ValT* A, int ADim1, int ADim2, const ValT* B, int BDim1, int BDim2, ValT* C)
{
if ( ADim2!=BDim1 ) throw std::runtime_error("Error sizes off");

memset((void*)C,0,sizeof(ValT)*ADim1*BDim2);
int cc2,cc1,cr1, a1,a2,a3;
for ( cc2=0 ; cc2<BDim2 ; ++cc2 ) {
    a1 = cc2*ADim2;
    a3 = cc2*BDim1
    for ( cc1=0 ; cc1<ADim2 ; ++cc1 ) {
          a2=cc1*ADim1;
          ValT b = B[a3+cc1];
          for ( cr1=0 ; cr1<ADim1 ; ++cr1 ) {
                    C[a1+cr1] += A[a2+cr1]*b;
           }
     }
  }
} 

Spróbuj, na pewno coś uratujesz.

Jeśli chodzi o twoje pytanie nr 1, powodem jest to, że mnożenie macierzy skaluje się jako O (n ^ 3), jeśli używasz trywialnego algorytmu. Istnieją algorytmy, które skalują się znacznie lepiej .

Stefano Borini
źródło
36
Przepraszam, ta odpowiedź jest całkowicie błędna. Implementacje BLAS nie są napisane w języku Fortran. Kod krytyczny dla wydajności jest napisany w asemblerze, a najbardziej popularne obecnie są napisane w C powyżej. BLAS określa również kolejność wierszy / kolumn jako część interfejsu, a implementacje mogą obsługiwać dowolną kombinację.
Andrew Tomazos
10
Tak, ta odpowiedź jest całkowicie błędna. Niestety jest pełen zdrowego bezsensu, np. Twierdzenie BLAS było szybsze dzięki Fortranowi. Posiadanie 20 (!) Pozytywnych ocen to zła rzecz. Teraz ten brak sensu rozprzestrzenia się nawet dalej ze względu na popularność Stackoverflow!
Michael Lehn
12
Myślę, że mylisz niezoptymalizowaną implementację referencyjną z wdrożeniami produkcyjnymi. Implementacja referencyjna służy jedynie do określenia interfejsu i zachowania biblioteki i została napisana w języku Fortran ze względów historycznych. Nie jest do użytku produkcyjnego. W produkcji ludzie używają zoptymalizowanych implementacji, które zachowują się tak samo jak implementacja referencyjna. Przestudiowałem wnętrze ATLAS (który wspiera Octave - Linux "MATLAB") i mogę potwierdzić, że jest wewnętrznie napisany w C / ASM. Niemal na pewno również komercyjne wdrożenia.
Andrew Tomazos
5
@KyleKanos: Tak, oto źródło ATLAS: sourceforge.net/projects/math-atlas/files/Stable/3.10.1 O ile wiem, jest to najczęściej używana przenośna implementacja BLASa typu open source. Jest napisany w C / ASM. Producenci wysokowydajnych procesorów, tacy jak Intel, również zapewniają implementacje BLAS, specjalnie zoptymalizowane dla ich układów. Gwarantuję, że na niskim poziomie części biblioteki Intels są napisane w (duuh) zestawie x86 i jestem prawie pewien, że części średniego poziomu byłyby napisane w C lub C ++.
Andrew Tomazos,
9
@KyleKanos: Jesteś zdezorientowany. Netlib BLAS jest implementacją referencyjną. Implementacja referencyjna jest znacznie wolniejsza niż implementacje zoptymalizowane (zobacz porównanie wydajności ). Kiedy ktoś mówi, że używa netlib BLAS w klastrze, nie oznacza to, że faktycznie używa implementacji referencyjnej netlib. To byłoby po prostu głupie. Oznacza to po prostu, że używają biblioteki z tym samym interfejsem co netlib blas.
Andrew Tomazos