C++Guns – RoboBlog blogging the bot

Die Performance ist ein Opfer von Level of Indirection. Jetzt auch bei mir. Ich will keinen langweiligen Code zeigen wie so etwas aussieht, oder wie man es verhindern kann. Viel interessanter ist es doch, wie das aus Assembler Ebene aussieht.
Dazu zwei Assembler Auszüge. Der erste mit zwei zusätzlichen Indirektionen. Es wird auf das i-th Element eines std::vector zugegriffen und zurück gegeben.

func(std::vector<T> const&, int):
        movq    %rdi, %rax       }
        movslq  %edx, %rdx       |
        salq    $4, %rdx         | Addressberechnung
        addq    (%rsi), %rdx     }
        movq    (%rdx), %rdx     * Indirektion 
        movq    (%rdx), %rdx     * Indirektion
        movdqu  (%rdx), %xmm0    }
        movups  %xmm0, (%rdi)    |
        movdqu  16(%rdx), %xmm1  | Kopieren der 4 double Variablen
        movups  %xmm1, 16(%rdi)  }
        ret

Ohne zweimalige Indirektion:

func2(std::vector<T> const&, int):
        movq    %rdi, %rax        }
        movslq  %edx, %rdx        |
        salq    $5, %rdx          | Addressberechnung
        addq    (%rsi), %rdx      }
        movdqu  (%rdx), %xmm0     }
        movups  %xmm0, (%rdi)     |
        movdqu  16(%rdx), %xmm1   | Kopieren der 4 double Variablen
        movups  %xmm1, 16(%rdi)   }
        ret

Es sind zwei Dinge zu erkennen.
Erstens gibt es zwei zusätzliche mov Anweisungen, welche eine Dereferenzierung eines Pointers bewirken. Die zusätzliche Arbeit den Befehl zu bearbeiten ist nicht das Problem. Vielmehr, dass die Daten nicht mehr hintereinander im RAM stehen. Es wird im RAM umher gesprungen. Dies geht bekanntlich auf die Performance, Cache misses u.s.w.

Bei dem konkreten Problem, aus dem dieses Beispiel hervorgegangen ist, lief das Programm nach dem Entfernen der zwei Indirektionen 10 mal schneller!

Zweitens ist zu erkennen, dass die zurück gegebenen Variablen aus 4 double Werten besteht, aber nur zwei Kopieraktionen zu erkennen sind. Von dem Register rdx, nach xmm0, nach rdi. Und noch einmal mit einem Offset von 16Byte (2 double Werte). Der Compiler hat erkannt, dass SIMD angewendet werden kann. Das Kopieren von zwei Werten mit einer Instruktion. Und das ohne zusätzlichen Bemühungen meinerseits C++ Code extra SIMD tauglich zu schreiben!

std::enable_if is dead; long live concepts!

Man beachte, dass die Requires clause NACH dem Funktionskopf kommen kann.

#include <type_traits>

template<typename T>
void func(T) requires std::is_integral_v<T> {
}

auto func2() {
    // error: cannot call function 'void func(T) requires  is_integral_v<T> [with T = double]'
    // constraints not satisfied  'is_integral_v<T>' evaluated to false
   func(2.0); // Nope
}

Mit -Ofast wird bei folgenden Code die Power Funktion weg optimiert und durch die dritte Wurzel ersetzt.

auto func(double MANNING2, double qNorm, double H ) {
    return MANNING2 * qNorm / (pow(H,(7.0/3.0)));
}

Im Assembler Code ist nur noch cbrt zu erkennen, kein pow mehr.

func(double, double, double):
        subq    $40, %rsp
        movsd   %xmm0, 8(%rsp)    # MANNING2
        movapd  %xmm2, %xmm0      # H nach xmm0 für cbrt
        movsd   %xmm1, 24(%rsp)   # qNorm
        movsd   %xmm2, 16(%rsp)   # H
        call    cbrt              # cbrt(H) nach xmm0
        movsd   16(%rsp), %xmm2   # xmm2 = H
        movsd   24(%rsp), %xmm1   # xmm1 = qNorm
        mulsd   8(%rsp), %xmm1    # xmm1 = MANNING2*qNorm
        addq    $40, %rsp 
        mulsd   %xmm2, %xmm2      # xmm2 = H*H
        mulsd   %xmm0, %xmm2      # xmm2 = H*H*cbrt(H)
        divsd   %xmm2, %xmm1      # xmm1 = MANNING2*qNorm /(H*H*cbrt(H))
        movapd  %xmm1, %xmm0
        ret

Der Assembler Code zurück nach C++ übersetzt gibt folgenden:

auto func(double MANNING2, double qNorm, double H ) {
    return MANNING2 * qNorm / (H*H*cbrt(H));
}

Die pow Funktion mit zwar bekannten, aber "unschönen" Exponenten wurde ersetzt durch eine Multiplikation und eine Power Funktion mit bekannte und "schönen" Exponenten: 1/3

Aber es könnte besser sein. Es fällt auf, dass die Funktionsargumente erst einmal auf den Stack gepusht werden. Offenbar braucht die cbrt Funktion die Register xmm0, xmm1, xmm2 selbst. Und so müssen die Werte eben gesichert werden. Da sie erst danach gebraucht werden. ...
Natürlich kann man MANNING2*qNorm berechnen bevor vor Nenner der Formel fertig berechnet ist. Mathematisch kann man das tun, die Klammern kann man setzten.

Aber den Compiler stört das nicht. Es wird immer erst cbrt vor MANNING2*qNorm. Ich denke cbrt braucht einfach mehr Zeit als eine Multiplikation. Und daher wird dieser Befehl als erstes ausgeführt. Und zwischenzeitlich können davon unabhängige Variablen bearbeitet werden.

Aber eine Sache kann dennoch optimiert werden. Indem H als ersten Parameter übergeben wird. Dann liegt es gleich im xmm0 Register für cbrt bereit und muss nicht erst dorthin kopiert werden. An der Anzahl der Instruktionen ändert dies nichts. Aber die mov Operation wird auch erst nach cbrt ausgeführt.

Sind Gleitkommazahlen sortierbar? Wegen NaN gibt es keine Ordnung und das Sortieren mit NaN im Datensatz schlägt fehlt.
Aber man kann sich schnell selbst eine Ordnung bauen:

#include <vector>
#include <algorithm>
#include <iostream>
#include <limits>
#include <cmath>

using namespace std;
using nl = numeric_limits<double>;

int main() {
  cout << boolalpha;

  vector<double> v { 1.0, 2.0, 0.0, 1.0/3.0, -0.0,
                    nl::quiet_NaN(), nl::infinity(), nl::epsilon(), nl::denorm_min(),
                    nl::max(), nl::min(), nl::lowest()
  };

  sort(v.begin(), v.end());

  for(const auto x : v) {
    cout << x << " ";
  }
  cout << "\nis sorted " << is_sorted(v.begin(), v.end()) << " << offensichtlich falsch\n";

  auto comp = [](const auto lhs, const auto rhs) {
    if(isnan(lhs)) return false;
    if(isnan(rhs)) return true;
    return lhs<rhs;
  };

  sort(v.begin(), v.end(), comp);

  for(const auto x : v) {
    cout << x << " ";
  }
  cout << "\nis sorted " << is_sorted(v.begin(), v.end()) << "\n";
}

-1.79769e+308 0 -0 0.333333 1 2 nan 4.94066e-324 2.22507e-308 2.22045e-16 1.79769e+308 inf 
is sorted true << offensichtlich falsch
-1.79769e+308 0 -0 4.94066e-324 2.22507e-308 2.22045e-16 0.333333 1 2 1.79769e+308 inf nan 
is sorted true << richtig

Why has forward_list no size() method?

In short: Its about performance. A single linked list size() function has a asymptotic complexity of O(n). We do not usually need something like that. So we don't have something.
The second reason is: For a O(1) size() function, it would increase the size of std::forward_list by 8 byte. We don't pay for what we don't need.

N2543 is the proposal, and it has a detailed discussion about size() [1]

The choice between Option 3 [not providing size()] and Option 2 [providing a O(1) size()] is more a matter of judgment. I have chosen Option 3 for the same reason that I chose insert-after instead of insert-before: Option 3 is more consistent with the goal of zero overhead compared to a hand-written C-style linked list. Maintaining a count doubles the size of a forward_list object (one word for the list head and one for the count), and it slows down every operation that changes the number of nodes. In most cases this isn't a change in asymptotic complexity (the one change in asymptotic complexity is in one of the forms of splice), but it is nonzero overhead. It's a cost that all users would have to pay for, whether they need this feature or not, and, for users who care about maintaining a count, it's just as easy to maintain it outside the list, by incrementing the count with every insert and decrementing it with every erase, as it is to maintain the count within the list.

But there is a very simple solution for this problem: A non-member size() function for std::forward_list as they are introduced in C++17.

namespace std {
    template<typename T>
    size_t size(const std::forward_list<T>& list) {
        size_t i=0;
        for(const auto& x : list) {
            (void)x;
            i++;
        }
        return i;
    }
}

But there is one problem: The std::size() now take O(1) or O(n) time, depending on the passed container.

[1] https://stackoverflow.com/questions/31822494/c-stl-why-has-forward-list-no-size-method

The problem seems to be an interesting little exercise that John Bentley once proposed to me: Insert a sequence of random integers into a sorted sequence, then remove those elements one by one as determined by a random sequence of positions: Do you use a vector (a contiguously allocated sequence of elements) or a linked list?

http://www.stroustrup.com/bs_faq.html

* Generate N random integers and insert them into a sequence so that each is inserted in its proper position in the numerical order. 5 1 4 2 gives:
- 5
- 1 5
- 1 4 5
- 1 2 3 4
* Remove elements one at a time by picking a random position in the sequence and removing the element there. Positions 1 2 0 0 gives
- 1 2 4 5
- 1 4 5
- 1 4
- 4
* For which N is it better to use a linked list than a vector (or an array) to represent the sequence?

GoingNative 2012 Bjarne Stroustrup: C++11 Style

Challenge accepted. To be continued...

Das restricted keyword in C

restrict keyword...is basically a promise to the compiler that for the scope of the pointer, the target of the pointer will only be accessed through that pointer (and pointers copied from it).

In C++ geben wir keine Versprechen. Grundsätzlich nicht.
Das Beispiel von cppreference.com verdeutlich das Optimierungspotential ziemlich anschaulich. Gegeben ist folgende Funktion:

int foo(int &a, int &b) {
    a = 5;
    b = 6;
    return a + b;
}

foo(int&, int&):
        movl    $5, (%rdi)     # store 5 in a
        movl    $6, (%rsi)     # store 6 in b
        movl    (%rdi), %eax   # read back from a in case previous store modified it
        addl    $6, %eax       # add 6 to the value read from a
        ret

Welche Optimierungsmöglichkeiten gibt es hier? Natürlich, der Rückgabewert ist 11 (5+6), solange die Referenzen von a und b auf unterschiedliche Variablen zeigen. Sonst ist der Rückgabewert 12 (6+6). Mit dem restrict keyword kann man lügen und behaupten, dass die Referenzen von a und b nicht auf die selbe Variablen zeigen. Aber so etwas haben wir in C++ nicht nötigt; der Compiler weiss es besser als wir.
Je nachdem in welchen Kontext die Funktion foo aufgerufen wird, kann optimiert werden, oder auch nicht. Beispielsweise ist in Funktion func zur Compilezeit sichergestellt, dass die Referenzen auf unterschiedliche Variablen zeigen. Entsprechend wird der Rückgabewert von 11 voraus berechnet und im Code abgelegt.

auto func(std::vector<int>& arr, size_t i) {
    return foo(arr[i], arr[i+1]);
}

func(std::vector<int, std::allocator<int> >&, unsigned long):
        movq    (%rdi), %rax
        addq    $1, %rsi             # temp increment i
        movl    $5, -4(%rax,%rsi,4)  # store 5 in a
        movl    $6, (%rax,%rsi,4)    # store 6 in b
        movl    $11, %eax            # the result is 11, a compile-time constant
        ret

Zeigen die Referenzen hingegen auf die selbe Variable, wie in der Funktion func2, erkennt auch dies der Compiler und berechnet einen anderen Rückgabewert: 12. Auch in diesem Fall muss die Addition zur Laufzeit nicht ausgeführt werden.

auto func2(std::vector<int>& arr, size_t i) {
    return foo(arr[i], arr[i]);
}

func2(std::vector<int, std::allocator<int> >&, unsigned long):
        movq    (%rdi), %rax
        movl    $6, (%rax,%rsi,4)  # store 6 in a and b.
        movl    $12, %eax          # the result is 12, a compile-time constant
        ret

Fazit: Hinweise an den Compiler mit dem restrict keyword in C sind in C++ absolut überflüssig und führen nur zu undefinierten Verhalten des Programms und damit zu schwer findbaren Bugs. Da die Debugzeit 1/3 bis gerne 2/3 der gesamten Entwicklungszeit ausmachen (in der Praxis, nicht in der Theorie), ist der Einsatz von C statt C++ unwirtschaftlich.

[1] https://en.cppreference.com/w/c/language/restrict
[2] https://stackoverflow.com/questions/776283/what-does-the-restrict-keyword-mean-in-c

Einfach N Zufallszahlen ziehen ist einfach. Eigentlich gibt es darüber nicht viel zu sagen, aber es ist ein schönes Beispiel für Semantik. (Dank an Ben für die Inspiration).

Wir müssen immer den Spagat zwischen Lesbarkeit und Performance machen. Manchmal sind die Verantwortlichkeiten der Funktionen nicht klar zu definieren. Soll die Funktion, welche N Zufallszahlen zieht, auch die Container erstellen in denen sie gespeichert werden? Das wäre für die Lesbarkeit förderlich, aber für die Performance schlecht. Da bei jedem Aufruf ein neuer Container initialisiert werden muss. Und wenn die Funktion einen bereits initialisierten Container übergeben bekommt und auch auf diesen Weg die Daten zurück gegeben werden, wäre der Lesefluss gestört. Vor allem muss beim ersten Aufruf von Hand ein neuer Container initialisiert werden.

Für beide Varianten gilt: Keine von Hand geschriebenen Schleifen. Sie enthalten einfach unnötige potentielle Fehler. Ein einfacher std::for_each Aufruf ist besser.
Hier mein Versuch beide Welten abzudecken. Durch die fehlenden Concepts in C++17 bläht sich der Template Code etwas auf. Aber dieses Problem löst sich ja von alleine.

// 7 loc to just call for_each... semantics is important
template<typename Container, typename PRNG, typename Distribution>
auto drawSamples(Container& data, PRNG& r, Distribution& U) {
    static_assert(std::is_arithmetic_v<typename Container::value_type>);
    static_assert(std::is_invocable_v<PRNG>, "PRNG is not invocable");
    static_assert(std::is_invocable_v<Distribution, PRNG&>, "Distribution is not invocable with PRNG");
    static_assert(std::is_same_v<typename Container::value_type, typename Distribution::result_type >, "Return type of the distribution and value type of container are not the same");

    std::for_each(std::begin(data), std::end(data), [&](auto& x){x = U(r);} );
    return data; // good idea? yes. consistent with the other drawSamples overload
                 // NO! bad feeling. pass by reference and return by value - it must be a copy somewhere. even with return optimazion.
}

template<typename PRNG, typename Distribution>
auto drawSamples(int N, PRNG& r, Distribution& U) {
    static_assert(std::is_invocable_v<PRNG>, "PRNG is not invocable");
    static_assert(std::is_invocable_v<Distribution, PRNG&>, "Distribution is not invocable with PRNG");    

    std::vector<typename Distribution::result_type> data(N);
    drawSamples(data, r, U);
    return data;
}


auto func() {
    std::minstd_rand r;
    std::uniform_int_distribution<int> U(1, 6);
    std::vector data = drawSamples(10, r, U);
    data = drawSamples(data, r, U);
    return data;
}

Bin gespannt ob sich dies in der Praxis bewährt.

update:
Ich glaube pass by reference und return by value ist immer eine Kopie. Und damit ineffizient. Selbst mit return Optimierung. Bin mir aber nicht sicher.
Vielleicht sieht das ganze mit Iteratoren und Ranges besser aus.

Bin Heute über diesen GCC Bugfix gestolpert GCC 8 stopped using RMW (Read Modify Write) instructions on x86. Es gibt also Assembler Instruktionen die in einem Zug einen Wert laden, den Wert ändern, und ihn wieder zurück schreiben. Dabei spielt es keine Rolle, ob dabei erst noch ein Pointer dereferenziert werden muss, oder die Variable per Referenz an die Funktion übergeben wird.

Im Link ist auch ein kleines Beispiel enthalten, welches ich noch weiter vereinfacht habe:

void f(void);

void h(int& x) {
    if(!--x) 
}      

Mit g++ -O1 -fpeephole2 (alle außer Version 8.1, welcher ja den BUG enthält) enthält man folgenden optimierten Assembler Code:

h(int&):
        subl    $1, (%rdi)
        je      .L14
        ret
.L14:
        call    f()
        ret

Wichtig ist die erste Zeile. Hier wird 1 minus der Wert auf welchen die Referenz zeigt gerechnet und das Ergebnis zurück geschrieben. Das ist also eine RMW Instruktion. Ohne peephole Optimierung sieht der erzeugte Code wie folgt aus:

h(int&):
        movl    (%rdi), %eax
        subl    $1, %eax
        movl    %eax, (%rdi)
        je      .L14
        ret
.L14:
        call    f()
        ret

Hier werden die Schritte einzeln mit je einer Assembler Anweisung abgearbeitet. Es gibt wohl etliche dieser peephole Optimierungen. Bei denen wird immer nur ein kurzer Codeabschnitt analysiert und durch einen effizienteren ersetzt.

P0184R0 Differing begin and end types in range-based for.

Sentinel - value to terminate a loop after it was read

Let me try something. A list implemented as array - Multiple lists stored in a single array.
This can be helpful if the lists are very short e.g. less than 10 items.
For this example, every list item is a single ID. Iterate over a list results in jumping around on the array.
I don't implement the part which insert new items in the list. Iterate over the lists is much more interesting.

Example:

#include <vector>
#include <iostream>

using namespace std;

const int noIDvalue = -1;

class MultiArrayList {
public:
    MultiArrayList() {
        // fill with some test data
        data.resize(10);
        data[0] = 9;
        data[1] = 3;
        data[2] = 1;
        data[3] = 8;
        data[4] = 7;
        data[5] = noIDvalue;
        data[6] = noIDvalue;
        data[7] = noIDvalue;
        data[8] = 5;
        data[9] = 4;

        startIDs.resize(3);
        startIDs[0] = 0;
        startIDs[1] = 2;
        startIDs[2] = 6;
    }


    std::vector<int> data;
    std::vector<int> startIDs;
};

int main() {
    MultiArrayList multiList;
    // print list content
    for(int nextID : multiList.startIDs) {
        cout << "List content: ";
        while(nextID != sentinelID) {
            cout << nextID << " ";
            nextID = multiList.data[nextID];
        }
        cout << "\n";
    }
}

List content: 0 9 4 7
List content: 2 1 3 8 5
List content: 6

The hand crafted while loop is very ugly. I want something with a range. e.g.

    for(size_t i=0; i < multiList.size(); ++i) {
        cout << "List content: ";
        for(const int ID : multiList.at(i)) {
            std::cout << ID << " ";
        }
        cout << "\n";
    }

Lets see if I can implement this with the new C++17.
We need MultiArrayList::at which simply returns a MultiArrayListRange. The range provide function begin() and end(). And we need MultiArrayListIterator which implements the iterator logic and of course MultiArrayListSentinel .

    auto MultiArrayList::at(size_t id) const {
        return MultiArrayListRange(startIDs.at(i), data);
    }

class MultiArrayListIterator {
public:
    MultiArrayListIterator(const int startID, const std::vector<int>& data)
        : curID(startID), data(data)
    {

    }

    void operator++() {
        curID = data[curID];
    }

    int operator*() const {
        return curID;
    }

private:
    int curID;
    const std::vector<int>& data;
};

template<int Delim>
struct MultiArrayListSentinel {

};

template<int Delim>
bool operator!=(const MultiArrayListIterator& lhs, const MultiArrayListSentinel<Delim>) {
    return *lhs != Delim;
}

class MultiArrayListRange {
public:
    MultiArrayListRange(const int startID, const std::vector<int>& data)
        : startID(startID), data(data)
    {

    }

    auto begin() {
        return MultiArrayListIterator(startID, data);
    }

    auto end() {
        return MultiArrayListSentinel<noIDvalue>();
    }

private:
    const int startID;
    const std::vector<int>& data;
};

Here we go. The core idea is the different return types of begin() and end(). One returns an iterator, the other a sentinel. The new operator!=() take both and check if the current ID is a noIDvalue, e.g. -1. If so, we are at the end of the list.

The end iterator is never incremented, decremented, or dereferenced. Requiring it to be an iterator serves no practical purpose.

As always: storing references is a bad idea. But everybody knows, iterators and now ranges are short-lived temporary objects which gets invalidated, so this is okay.

PS: operator!= with sentinel is faster, because the right hand side is a compile time constant. It saves one assembler instruction (25%).