Kubin

• Strona główna • Moje rzeczy • Kartki • Esoterica • Kontakt •

Tajemnice libstdc++

Disclaimer

Artykuł oczywiście dotyczy stricte C++. Moja wersja g++ w chwili pisania kartki to 5.3.0. Przed wyciąganiem wniosków warto sprawdzić własną wersję kompilatora lub ustalenia techniczne sprawdzarki. W każdym razie, do swojej wersji będę się tutaj odnosił. Będę też zakładał że korzystamy ze standardu C++11. Dodatkowo będę wycinał komentarze z wklejanego kodu źródłowego.
Ogólnie cała ta kartka jest pełna technikaliów. Casuals beware.

libstdc++?

Każda wersja g++ ma w pakiecie implementację biblioteki C++, czyli STL. Jednak to nie ten sam zespół jest odpowiedzialny za tę implementację, co za kompilator - owa jest oddzielnym projektem o nazwie libstdc++. Interesują nas pliki źródłowe, których trzeba poszukać w swojej instalacji. W moim MinGW (które pewnie będzie trochę zbliżone do innych instalacji g++) są w /lib/gcc/mingw32/5.3.0/include/c++. Względem takiej ścieżki, w której znajdują się wszystkie pliki znacznikowe (algorithm i tym podobne), będę się tutaj odnosił.

Wprowadzenie

Przede wszystkim, po chwili poszukiwań, znajdziemy tu implementację każdej rzeczy z STLa, taką jaką obiecuje dokumentacja. Należy wziąć pod uwagę, że kod, który ujrzymy, będzie bardzo mocno opatrzony szablonami, makrami związanymi z implementacją kompilatora, makrami od debugu, etc., przez co może być trochę nieczytelny. Tak na przykład wygląda implementacja std::reverse (algorithm odeśle nas do bits/stl_algo.h oraz bits/stl_algobase.h, przy czym szukana implementacja jest w tym pierwszym):

std::reverse 
template<typename _BidirectionalIterator>
    inline void
    reverse(_BidirectionalIterator __first, _BidirectionalIterator __last)
    {
      // concept requirements
      __glibcxx_function_requires(_Mutable_BidirectionalIteratorConcept<
				  _BidirectionalIterator>)
      __glibcxx_requires_valid_range(__first, __last);
      std::__reverse(__first, __last, std::__iterator_category(__first));
    }

Widać tutaj to, o czym mówiłem. Szablon – ✓. Makra – ✓. Debug: – ✓. No i do tego odwołanie do __reverse, aby ukryć implementację (podwójny _ często oznacza, że coś jest używane przez wewnętrzną implementację i lepiej tego nie psuć, podobna konwencja dotyczy pojedynczego). Patrzymy na tę funkcję i widzimy:

__reverse 
template<typename _BidirectionalIterator>
    void
    __reverse(_BidirectionalIterator __first, _BidirectionalIterator __last,
	      bidirectional_iterator_tag)
    {
      while (true)
	if (__first == __last || __first == --__last)
	  return;
	else
	  {
	    std::iter_swap(__first, __last);
	    ++__first;
	  }
    }

    template<typename _RandomAccessIterator>
    void
    __reverse(_RandomAccessIterator __first, _RandomAccessIterator __last,
	      random_access_iterator_tag)
    {
      if (__first == __last)
	return;
      --__last;
      while (__first < __last)
	{
	  std::iter_swap(__first, __last);
	  ++__first;
	  --__last;
	}
    }

Cóż, to prawie to czego szukamy. Pozostaje rozszyfrować czym jest std::iter_swap. Ja nie wiedziałem. Cóż, C++ nie przestaje zaskakiwać. Tym razem znajdziemy go w bits/stl_algobase.h

std::iter_swap 
template<bool _BoolType>
    struct __iter_swap
    {
      template<typename _ForwardIterator1, typename _ForwardIterator2>
        static void
        iter_swap(_ForwardIterator1 __a, _ForwardIterator2 __b)
        {
          typedef typename iterator_traits<_ForwardIterator1>::value_type
            _ValueType1;
          _ValueType1 __tmp = _GLIBCXX_MOVE(*__a);
          *__a = _GLIBCXX_MOVE(*__b);
          *__b = _GLIBCXX_MOVE(__tmp);
	}
    };

  template<>
    struct __iter_swap<true>
    {
      template<typename _ForwardIterator1, typename _ForwardIterator2>
        static void 
        iter_swap(_ForwardIterator1 __a, _ForwardIterator2 __b)
        {
          swap(*__a, *__b);
        }
    };

Pozwoliłem sobie oszczędzić pokazywania 50 linii boilerplate'u, który faktycznie nazywa się std::iter_swap (A wy napisaliście kiedyś odwracanie vectora w 150 linii?). Chodzi o prosty wniosek: znalezienie czegoś tutaj może być... skomplikowane. Można jednak znaleźć kwiatki.

Rozszerzenia GNU

Owe rozszerzenia znajdują się głównie w ext, i głównie o będę opowiadał. Większość ze znajdujących się tam rzeczy jest w przestrzeni nazw __gnu_cxx. Czeka na nas tutaj trochę ciekawych rzeczy, o ile jesteśmy dostatecznie uparci.

Szybkie potęgowanie

W pliku ext/numeric znajdziemy __gnu_cxx::power - implementację szybkiego potęgowania działającą w \(O(\log{n})\). Wywołanie jest proste: power(a, b) zwróci a podniesione do potęgi b. Dodatkowo można przekazać trzeci parametr będący obiektem który można wywołać (tj. ze zdefiniowanym operator()). Ten obiekt (oznaczmy go F) musi też mieć zdefiniowany wynik funkcji identity_element(F), która posłuży za wynik dla power(a, 0). Wtedy dokładniejsza definicja power to:

Definicja power 
def power(a, b, F):
    if b == 0:
        return identity_element(F)
    else:
        return F(a, power(a, b-1, F))

Oczywiście power ma sens dla całkowitych i nieujemnych b.

Często będziemy chcieli liczyć potęgi modulo pewna liczba, i wtedy będziemy musieli napisać tzw. functor, czyli prosty typ, który definiuje operator(). Dużo przykładów takich prostych functorów znajdziemy w <functional>, domyślnym w tym wypadku jest std::multiplies<T>. Inne przykłady: std::plus<T>, std::modulus<T>.

Taki functor (wraz z konstruktorem do dowolnego modulo, warto jednak modulować przez stałą, patrz: O modulo) oraz przykład wykorzystania w __gnu_cxx::power, może wyglądać tak:

Przykład użycia __gnu_cxx::power 
#include <iostream>
#include <ext/numeric>

using namespace __gnu_cxx;
using namespace std;

template<typename T>
struct mod_multiplies
{
    T mod;
    mod_multiplies(T _mod) : mod(_mod) {}
    T operator() (T a, T b) { return (a*b) % mod; }
};
template<typename T>
T identity_element(mod_multiplies<T>) { return 1; }

int main()
{
    int64_t a, m; uint64_t b;
    cin >> a >> b >> m;
    cout << power(a, b, mod_multiplies<int64_t>(m));
}

Jeżeli ktoś bardzo chce skrócić kod, można: 1) nie robić szablonu, 2) odziedziczyć od multiplies, aby wartość identity_element była zdefiniowana, 3) usunąć konstruktor i modulować przez stałą.

Ciekawostka: kompilator umie zoptymalizować proste, stałe wywołania jak power(3, 5) do stałej (243), albo power(x, 3) do dosłownego x*x*x. Dobrze mieć to na uwadze, zamiast pisać to samo wiele razy.

Według mnie opłaca się znać – mniej miejsca na błąd w implementacji, a functor napisać bardzo łatwo jak się wyuczy. Jeżeli zdefniujemy dla naszego typu operator* nie trzeba będzie pisać dodatkowego functora.

Można w ten sposób potęgować macierze (trzeba tylko napisać odpowiedni identity_element – macierz tożsamościową) i inne operacje, które sensownie wpasowują się w definicję.

Policy-Based Data Structures

Notatka o notacji \(O\)

Notacja \( g(x) = \Theta(f(x)) \) oznacza, że pewna funkcja \(g(x)\) jest asymptomatycznie ograniczona przez \(f(x)\) zarówno od dołu jak i od góry. Samo \(O\) zwykle oznacza ograniczenie od góry na pewną funkcję – w informatyce zwykle maksymalny bądź średni czas działania algorytmu. Więcej do poczytania o notacji \(O\) w bibliografii.

Jeden z najbardziej rozbudowanych i najmniej udokumentowanych bytów jakie spotkałem, policy-based data structures są bardzo ciekawą rzeczą która może nas uratować od utopienia się w implementacji jakiejś struktury, takiej jak BST (Binary Search Tree [binarne drzewo poszukiwań] - nawiasem mówiąc temat ciekawy, który zasługuje na własną kartkę) czy hashmapa. Bez wątpienia - warto znać. Czemu policy-based? Ponieważ ich działanie możemy modyfikować wybierając tag (czyli jakiej dokładnie struktury danych użyć aby zaimplementować dany szablon) oraz inne fragmenty decydujące o działaniu struktury. Wszystkie są w przestrzeni nazw __gnu_pbds.

Dzielą się na:

Pisanie zmodyfikowanych BST za pomocą tree

Czasami implementując nasze rozwiązanie możemy dojść do wniosku, że zapewniona przez STL implementacja BST (std::map albo std::set) może być niewystarczająca (Ot, taki przykład, zadanie "Rotacje na drzewie" z XVIII OI). Wtedy trzeba by samemu napisać takie BST. Dwa konkretne przykłady bardziej skomplikowanych drzew, które mogą nam być potrzebne:

  • Order statistics tree – drzewo, które pozwala na zliczanie, ile elementów jest mniejszych od danej wartości, oraz na znajdowanie i-tego elementu w posortowanej kolejności.
  • Interval tree – drzewo, które pozwala na operacje na zbiorze przedziałów. Na przykład można pytać o ilość przedziałów przecinających dany punkt czy inny przedział.

    • Oczywiście czasami takie bardziej skomplikowane BST można symulować zwykłym drzewem przedziałowym czy Fenwickiem. Jednakże należy wziąć pod uwagę, że takie rozwiązanie może być nieoptymalne albo skomplikowane w implementacji. (w przytoczonym zadaniu należy samemu napisać merge na dynamicznych drzewach przedziałowych)

    Ale! Możemy tak zmodyfikować tree, aby zapewniało nam potrzebną funkcjonalność. Przedstawię, jak to zrobić na przykładzie order statistics tree.

    Aby implementować nasze rozszerzenia wykorzystamy dostarczony przez tree interfejs NodeUpdate.

    Order statistics tree

    Przede wszystkim, order statistics tree jest już zaimplementowane za pomocą (jedynego) wbudowanego NodeUpdate, czyli tree_order_statistics_node_update. Przykład użycia:

    #include <bits/stdc++.h>
#include <ext/pb_ds/assoc_container.hpp>

using namespace std;
using namespace __gnu_pbds;

template<typename T>
using indexed_set = tree<
    T, null_type, less<T>, rb_tree_tag, tree_order_statistics_node_update
>;

int main()
{
    indexed_set<int> S;
    S.insert(1); S.insert(4); S.insert(7);
    cout << *S.find_by_order(1) << endl; // 4
    cout << S.order_of_key(0) << " " << S.order_of_key(4) << endl; // 0 1
}

    Jednak jak wziąć się za implementację własnej aktualizacji wierzchołków? Należy przeanalizować tę wbudowaną. Idea NodeUpdate opiera się na przechowywaniu metadanych w wierzchołkach, tak, że metadane można wywnioskować tylko na podstawie wierzchołka oraz jego dzieci (niestety z powodu tego ograniczenia nie można pisać np. własnej opóźnionej propagacji, co można robić za pomocą napisanego własnoręcznie BST). Aby udostępniać jakąś funkcjonalność, musimy też napisać odpowiednie funkcje. Żeby móc tego dokonać, dostajemy dostęp do wewnętrznego interfejsu drzewa.

    Najpierw zastanówmy się jak będzie wyglądać implementacja. Metadanymi będzie ilość wierzchołków w poddrzewie. Wtedy ilość wierzchołków w lewym poddrzewie to ilość elementów mniejszych. Na podstawie tego można wymyślić wszystkie potrzebne nam operacje.

    A oto przykład implementacji order statistics tree (używa się jej tak jak wbudowanej):

    Implementacja order statistics tree 
    #include <bits/stdc++.h>
#include <ext/pb_ds/assoc_container.hpp>

using namespace std;
using namespace __gnu_pbds;
using __gnu_pbds::detail::branch_policy;

template<typename node_const_iterator, typename node_iterator,
         typename cmp_fn, typename _Alloc>
struct custom_tree_order_statistics_node_update
    : branch_policy<node_const_iterator, node_iterator, _Alloc>
{
    typedef size_t metadata_type;

    typedef typename node_iterator::value_type iterator;
    typedef typename node_const_iterator::value_type const_iterator;
    typedef typename iterator_traits<iterator>::value_type key_type;

    virtual node_iterator node_begin() = 0;
    virtual node_iterator node_end() = 0;
    virtual node_const_iterator node_begin() const = 0;
    virtual node_const_iterator node_end() const = 0;
    virtual cmp_fn& get_cmp_fn() = 0;
    virtual ~custom_tree_order_statistics_node_update() {}

    inline void operator()(node_iterator it, node_const_iterator end_it) const
    {
        auto L = it.get_l_child(); auto R = it.get_r_child();
        const_cast<metadata_type&>(it.get_metadata()) =
            (L == end_it ? 0 : const_cast<metadata_type&>(L.get_metadata())) +
            (R == end_it ? 0 : const_cast<metadata_type&>(R.get_metadata())) + 1;
    }
    iterator find_by_order(size_t order)
    {
        auto it = node_begin(); auto end_it = node_end();
        while(it != end_it)
        {
            auto L = it.get_l_child();
            size_t oL = L == end_it ? 0 : L.get_metadata();
            if(oL == order)
                break;
            else if(oL > order)
                it = L;
            else
            {
                it = it.get_r_child();
                order -= oL + 1;
            }
        }
        return *it;
    }
    size_t order_of_key(const key_type& key)
    {
        auto it = node_begin(); auto end_it = node_end();
        const cmp_fn& cmp = get_cmp_fn();
        size_t order = 0;
        while(it != end_it)
        {
            node_const_iterator L = it.get_l_child();
            if(cmp(key, this->extract_key(**it)))
                it = L;
            else if(cmp(this->extract_key(**it), key))
            {
                order += 1 + ((L == end_it) ? 0 : L.get_metadata());
                it = it.get_r_child();
            }
            else
            {
                order += (L == end_it) ? 0 : L.get_metadata();
                break;
            }
        }
        return order;
    }
};
template<typename T, typename Tm = null_type>
using order_statistics_tree = tree<
    T, Tm, less<T>, rb_tree_tag, custom_tree_order_statistics_node_update
>;

    Rozszerzenia kompilatora

    Notatka: większość funkcji postaci __builtin_[...] jest dostępna globalnie i nie wymaga żadnego pliku nagłówkowego.

    Mimo tego, że rozszerzenia kompilatora nie mają zbyt dużo wspólnego z samym libstdc++, uznałem, że tematycznie tutaj pasują.

    Zaawansowane operacje bitowe

    Czasami potrzebujemy czegoś bardziej skomplikowanego niż zwykłe operacje bitowe. W GCC do dyspozycji dostajemy bliźniacze funkcje: zliczanie zer wiodących oraz kończących w reprezentacji bitowej liczby. Funkcje te (kryją się pod nazwami int __builtin_clz(int x) (count leading zeros) oraz int __builtin_ctz(int x) (count trailing zeros)). Należy wziąć pod uwagę, że zachowanie funkcji dla 0 jest niezdefiniowane. Wynika to z tego, że zwykle te funkcje kompilator zapisze jako wywołanie odpowiedniej instrukcji procesora (na systemach gdzie ta instrukcja nie jest dostępna – kompilator ma swoją zastępczą implementację), a różne procesory mogą różnie obsługiwać ten przypadek.

    Ponieważ domyślnie działają na int, funkcje występują w odmianach pracujących na long i long long, i aby je wydobyć należy (w stylu C) do nazwy funkcji dodać sufiks l lub ll.

    Do czego mogą nam się przydać te operacje? Przede wszystkim clz okaże się bardzo użyteczne. Okazuje się, że dla danej dodatniej liczby \( n \) indeks najbardziej znączącej jedynki (indeks, czyli wykładnik potęgi, którą reprezentuje) jest jednocześnie równy \( \lfloor \log_2{n} \rfloor \). A skąd wziąć ten indeks? Jest to oczywiście ilość wszystkich bitów (minus jeden) odjąć ilość zer wiodących. Na podstawie tego możemy sobie napisać taką oto funkcję:

    uint32_t log2floor(uint32_t x)
{ 
    return 31 - __builtin_clz(x); 
}

    A jak otrzymać sufit z logarytmu dwójkowego, który też czasami jest potrzebny? (Na przykład przy drzewie przedziałowym - w ten sposób otrzymamy wysokość tudzież ilość liści). Tym razem okazuje się, że dla liczb całkowitych spełniona jest tożsamość \( \lceil \log_2{n} \rceil = \lfloor \log_2{(2n-1)} \rfloor \) (Można to sobie uzasadnić techniką 'dodaj jak najwięcej, ale nie tak, żeby zmienić wynik', której używamy przy dzieleniu z sufitem zamiast podłogą - \( \lceil \frac{a}{b} \rceil = \lfloor \frac{a+b-1}{b} \rfloor \)).

    Słowem komentarza: funkcja log2floor jest dostępna jako std::__lg w <algorithm>. Jest ona szablonem poprawnie obsługującym wszystkie typy całkowite co najmniej wielkości int – warto pamiętać.

    Poza zliczaniem zer można też zliczać wszystkie jedynki. Funkcja ta, zwana za morzem popcount jest nam dana jako int __builtin_popcount(int x). Także działa na int i ma odmiany do innych typów z odpowiednimi sufiksami. Analogicznie też jest zastąpywana instrukcją, o ile jest to możliwe. Przydaje się np. do liczenia wielkości pewnego podzbioru przy pracy na bitmaskach.

    Wielkie nadzieje

    Poza w miarę logicznymi operacjami bitowymi, mamy też specjalne funkcje dające kompilatorowi więcej informacji, które mają mu ułatwić optymalizację kodu. Z tych łatwiejszych do użycia mamy long __builtin_expect(long value, long expected) (wraz z long __builtin_expect_with_probability(long value, long expected, double probability)) oraz void __builtin_unreachable().

    Unreachable wywołujemy, gdy chcemy kompilator zapewnić, że program nigdy nie dojdzie do pewnego miejsca w kodzie. Na przykład, jeżeli chcemy znaleźć w pewnym ciągu liczbę 0, a jesteśmy pewni że owa w nim występuje, i przechodzimy pętlą, returnując indeks jeżeli znaleźliśmy owe 0, to po pętli można dać __builtin_unreachable(); - bo zawsze returnujemy przed zakończeniem pętli.

    Natomiast expect wykorzystujemy, aby dać kompilatorowi wskazówkę, jak zoptymalizować pewien wyjątek prowadzący do rozgałęzienia (tak zwany branching). Na przykład, jeżeli jesteśmy pewni, że warunek zachodzi właściwie zawsze (a jednak czasami nie), to używamy właśnie __builtin_expect. Zwróci nam niezerową wartość, jeżeli warunek zachodzi.
    Tematyczny przykład: w funkcji log2floor, aby zdefiniować wartość wywołania dla 0 można użyć __builtin_expect(x != 0, true) ? 31 - __builtin_clz(x) : -1u, zakładając, że funkcja prawie nigdy nie będzie wywoływana z 0. W porównaniu z x != 0 ? 31 - __builtin_clz(x) : -1u: kompilator wygeneruje kod, w którym procesor, zamiast domyślnie najpierw ładować na stos -1u, będzie preferował liczenie clz, oraz zamiast jumpować do miejsca, gdzie liczone jest clz, będzie jumpował do miejsca gdzie zwracane jest -1u. Teoretycznie rzecz biorąc, kod został lepiej zoptymalizowany - jumpy są kosztowne, lepiej żeby były wykorzystywane w przypadku rzadszym. Inny przykład to przy modulo przez liczby Mersenne'a ('O modulo'): można w ten sposób wyifować przypadek, kiedy \(x\) jest wielokrotnością \(2^k - 1\).
    Domyślne wywołanie zakłada, że warunek zachodzi w 90% przypadków. Jeżeli chcemy założyć inne prawdopodobieństwo, należy użyć __builtin_expect_with_probability (w analogiczny sposób jak wersję zwykłą), i podać prawdopodobieństwo jako liczbę od 0 do 1. Podane prawdopodobieństwo pozwala kompilatorowi na ocenę, jaką optymalizację zastosować.

    Bibliografia