Warum brauchen wir Bereiche von C ++ 20 in einem einfachen Brecher?

In letzter Zeit wurden die Bereiche, die in den C ++ 20-Standard aufgenommen werden sollten, viel diskutiert, einschließlich Habré ( ein Beispiel, in dem es viele Beispiele gibt ). Es gibt genug Kritiker von Intervallen, das sagen sie

• Sie sind zu abstrakt und werden nur für sehr abstrakten Code benötigt
• Die Lesbarkeit des Codes mit ihnen verschlechtert sich nur
• Intervalle verlangsamen den Code

Mal sehen, ganz Arbeiter-Bauer praktische Aufgabe, um zu verstehen, ob diese Kritik gültig ist und ist es wahr, dass Eric Nibler von Bartosz Milewski gebissen wurde und Range-v3 nur bei Vollmond schreibt .

Wir werden die folgende Funktion mit der Trapezmethode integrieren:$$ $ inline $ f (t) = 3 t ^ 2 \ sin t ^ 3 $ inline $ im Bereich von Null bis$$ $ inline $ \ tau $ inline $ . Wenn$$ $ inline $ \ tau ^ 3 / \ pi $ inline $ gleich einer ungeraden Zahl, dann ist das Integral 2.

Das Problem: Wir schreiben einen Prototyp einer Funktion, die das Integral nach der Trapezmethode berechnet . Auf den ersten Blick scheint es, dass hier keine Abstraktionen benötigt werden, aber Geschwindigkeit ist wichtig. In der Tat ist dies nicht ganz richtig. Für die Arbeit muss ich oft "Zahlenbrecher" schreiben, deren Hauptnutzer ich selbst bin. Also muss ich auch ihre Fehler unterstützen und damit umgehen (leider meine Kollegen - nicht immer nur ich). Und es kommt auch vor, dass Code nicht verwendet wird, beispielsweise ein Jahr, und dann ... Im Allgemeinen müssen auch Dokumentation und Tests geschrieben werden.

Welche Argumente sollte eine Integratorfunktion haben? Integrierbare Funktion und Gitter (Punktmenge$$ $ inline $ t_1, t_2, t_3 ... $ inline $ zur Berechnung des Integrals verwendet). Und wenn mit der integrierten Funktion alles klar ist ( std::function ist hier genau richtig), in welcher Form soll dann das Netz übertragen werden? Mal sehen.

Optionen

Zunächst schreiben wir - damit es etwas gibt, mit dem man die Leistung vergleichen kann - eine einfache for Schleife mit einem konstanten Zeitschritt:

 // trapezoidal rule of integration with fixed time step; // dt_fixed is the timestep; n_fixed is the number of steps double integrate() { double acc = 0; for(long long i = 1; i < n_fixed - 1; ++i) { double t = dt_fixed * static_cast<double>(i); acc += dt_fixed * f(t); } acc += 0.5 * dt_fixed * f(0); acc += 0.5 * dt_fixed * f(tau); return acc; } 

Wenn Sie diesen Zyklus verwenden, können Sie den Beginn und das Ende des Integrationsintervalls sowie die Anzahl der Punkte für diese Integration selbst als Argumente an die Funktion übergeben. Stop - die Trapezmethode findet auch mit einem variablen Schritt statt, und unsere integrierbare Funktion fordert einfach die Verwendung eines variablen Schritts an! Lassen Sie uns noch einen Parameter haben ($$ $ inline $ b $ inline $ ) um die "Nichtlinearität" zu kontrollieren und unsere Schritte zum Beispiel wie folgt aussehen zu lassen:$$ $ inline $ \ Delta t (t) = \ Delta t_0 + bt $ inline $ . Dieser Ansatz (Einführung eines zusätzlichen numerischen Parameters) wird wahrscheinlich an einer Million Stellen verwendet, obwohl sein Fehler anscheinend für alle offensichtlich sein sollte. Und wenn wir eine andere Funktion haben? Und wenn wir irgendwo in der Mitte unseres numerischen Intervalls einen kleinen Schritt brauchen? Was aber, wenn eine integrierbare Funktion mehrere Funktionen aufweist? Im Allgemeinen müssen wir in der Lage sein, jedes Gitter zu vermitteln. (Trotzdem werden wir in den Beispielen bis zum Ende die reale Trapezmethode „vergessen“ und der Einfachheit halber ihre Version mit einem konstanten Schritt betrachten, wobei zu berücksichtigen ist, dass das Gitter beliebig sein kann).

Da das Raster beliebig sein kann, übergeben wir seine Werte$$ $ inline $ t_1, t_2, ... $ inline $ eingewickelt in std::vector .

 // trapezoidal rule of integration with fixed time step double integrate(vector<double> t_nodes) { double acc = 0; for(auto t: t_nodes) { acc += dt_fixed * f(t); } acc -= 0.5 * dt_fixed * f(0); acc -= 0.5 * dt_fixed * f(tau); return acc; } 

Es gibt mehr als genug Communities in diesem Ansatz, aber was ist mit der Leistung? Mit Speicherverbrauch? Wenn früher alles auf dem Prozessor zusammengefasst wurde, müssen wir jetzt zuerst den Speicherbereich ausfüllen und dann daraus lesen. Und die Kommunikation mit dem Gedächtnis ist ziemlich langsam. Und die Erinnerung ist immer noch nicht Gummi ( und Silikon )

Schauen wir uns die Wurzel des Problems an. Was braucht ein Mensch, um glücklich zu sein? Was braucht unser Zyklus (bereichsbasiert für Schleife)? Jeder Container mit Iteratoren begin() und end() sowie ++ , * und != Operatoren für Iteratoren. Also werden wir schreiben.

 //   ,    ,      template <typename T> double integrate(T t_nodes) { double acc = 0; for(auto t: t_nodes) { acc += dt_fixed * f(t); } acc -= 0.5 * dt_fixed * f(0); acc -= 0.5 * dt_fixed * f(tau); return acc; } // ... //      class lazy_container { public: long long int n_nodes; lazy_container() { n_nodes = n_fixed; } ~lazy_container() {} class iterator { public: long long int i; // index of the current node iterator() { i = 0; } ~iterator() {} iterator& operator++() { i+= 1; return *this; } // ! bool operator!=(const iterator& rhs) const { return i != rhs.i; } double operator* () const { return dt_fixed * static_cast<double>(i); } }; iterator begin() { return iterator(); } iterator end() { iterator it; it.i = n_nodes; return it; } }; // ... //      lazy_container t_nodes; double res = integrate(t_nodes); 

Wir berechnen hier einen neuen Wert.$$ $ inline $ t_i $ inline $ auf Anfrage, genau wie wir es in einer einfachen for Schleife getan haben. Es gibt keine Speicherzugriffe und es besteht die Hoffnung, dass moderne Compiler den Code sehr effizient vereinfachen. Gleichzeitig hat sich der Code der Integrationsfunktion nicht wesentlich geändert und kann std::vector .

Wo ist die Flexibilität? Tatsächlich können wir jetzt jede Funktion in den ++ Operator schreiben. Das heißt, dieser Ansatz ermöglicht es tatsächlich, eine Funktion anstelle eines einzelnen numerischen Parameters zu übertragen. Das im laufenden Betrieb erzeugte Raster kann beliebig sein, und wir verlieren (wahrscheinlich) auch nicht an Leistung. Es fühlt sich jedoch überhaupt nicht lazy_container jedes Mal einen neuen lazy_container zu schreiben, um das Raster auf eine neue Weise zu verzerren (es sind die gleichen 27 Zeilen!). Natürlich können Sie die Funktion, die für die Erzeugung des Rasters verantwortlich ist, zu einem Parameter unserer Integrationsfunktion machen und lazy_container , lazy_container entschuldigen Sie, kapseln Sie ihn.

Sie fragen - dann wird wieder etwas falsch sein? Ja! Zunächst muss die Anzahl der zu integrierenden Punkte separat übertragen werden, was zu einem Fehler führen kann. Zweitens muss das erstellte nicht standardmäßige Fahrrad von jemandem unterstützt und möglicherweise entwickelt werden. Können Sie sich zum Beispiel sofort vorstellen, wie Sie mit diesem Ansatz einen Kombinator für Funktionen im ++ Operator ++ können?

C ++ seit über 30 Jahren. Viele in diesem Alter haben bereits Kinder, und C ++ hat nicht einmal Standard-Lazy-Container / Iteratoren. Ein Albtraum! Aber alles (im Sinne von Iteratoren, nicht von Kindern) wird sich bereits im nächsten Jahr ändern - der Standard (möglicherweise teilweise) wird die Range-v3- Bibliothek enthalten, die von Eric Nibler seit mehreren Jahren entwickelt wurde. Wir werden die Werke seiner Früchte verwenden. Der Code sagt alles für sich:

 #include <range/v3/view/iota.hpp> #include <range/v3/view/transform.hpp> //... auto t_nodes = ranges::v3::iota_view(0, n_fixed) | ranges::v3::views::transform( [](long long i){ return dt_fixed * static_cast<double>(i); } ); double res = integrate(t_nodes); 

Die Integrationsfunktion bleibt gleich. Das heißt, nur 3 Zeilen, um unser Problem zu lösen! Hier erzeugt iota_view(0, n) ein verzögertes Intervall (Bereich, ein Objekt, das den verallgemeinerten Anfang und das verallgemeinerte Ende kapselt; ein verzögerter Bereich ist eine Ansicht), das bei Iteration bei jedem Schritt die nächste Zahl im Bereich [0, n] berechnet. Es ist lustig, dass der Name ι (der griechische Buchstabe iota) vor 50 Jahren auf die APL-Sprache verweist . Stick | Mit dieser Option können Sie Pipelines von Intervallmodifikatoren schreiben, und die transform ist in der Tat ein solcher Modifikator, der mithilfe einer einfachen Lambda-Funktion eine Folge von Ganzzahlen in eine Reihe transform$$ $ inline $ t_1, t_2, ... $ inline $ . Alles ist einfach wie in ein Märchen Haskell.

Aber wie macht man einen variablen Schritt? Alles ist genauso einfach:

 #include <range/v3/view/drop.hpp> #include <range/v3/view/iota.hpp> #include <range/v3/view/transform.hpp> //... // trapezoidal rule of integration; step size is not fixed template <typename T> double integrate(T t_nodes) { double acc = 0; double t_prev = *(t_nodes.begin()); double f_prev = f(t_prev); for (auto t: t_nodes | ranges::v3::views::drop(1)) { double f_curr = f(t); acc += 0.5 * (t - t_prev) * (f_curr + f_prev); t_prev = t; f_prev = f_curr; } return acc; } //... auto step_f = [](long long i) { if (static_cast<double>(i) <= 1 / a) { return pow(2 * M_PI, 1/3.0) * a * static_cast<double>(i); } else { return tau * pow(static_cast<double>(i) / static_cast<double>(n), 1/3.0); } }; auto t_nodes = ranges::v3::iota_view(0, n) | ranges::v3::views::transform(step_f); double res = integrate(t_nodes); 

Ein aufmerksamer Leser stellte fest, dass wir in unserem Beispiel durch den variablen Schritt die Anzahl der Gitterpunkte um den Faktor drei reduzieren konnten, während zusätzlich spürbare Kosten für die Berechnung anfallen$$ $ inline $ t_i $ inline $ . Aber wenn wir noch einen nehmen$$ $ inline $ f (t) $ inline $ kann sich die Anzahl der Punkte viel mehr ändern ... (aber hier wird der Autor schon faul).

Also Timings

Wir haben folgende Möglichkeiten:

• v1 - einfache Schleife
• v2 -$$ $ inline $ t_i $ inline $ liegen in std::vector
• v3 - notdürftiger lazy_container mit notdürftigem Iterator
• v4 - Intervalle von C ++ 20 (Bereiche)
• v5 - Bereiche wieder, aber nur hier wird die Trapezmethode mit einer variablen Tonhöhe geschrieben

Hier ist, wofür es sich herausstellt (in Sekunden)$$ $ inline $ \ tau = (10 \, 000 \, 001 \ times \ pi) ^ {1/3} $ inline $ , für g ++ 8.3.0 und clang ++ 8.0.0 auf Intel® Xeon® CPU® X5550 (Anzahl der Schritte ungefähr$$ $ inline $ 1.5 \ times 10 ^ 8 $ inline $ , mit Ausnahme von v5, wo die Schritte dreimal kleiner sind (das Ergebnis der Berechnungen aller Methoden unterscheidet sich von den beiden um nicht mehr als 0,07):

Im Allgemeinen ging die Fliege über das Feld, die Fliege fand eine Münze!

Zusammenfassung

Selbst in einer sehr einfachen Aufgabe erwiesen sich Bereiche als sehr nützlich: Anstelle von Code mit selbst erstellten Iteratoren in mehr als 20 Zeilen haben wir 3 Zeilen geschrieben, ohne Probleme mit der Lesbarkeit des Codes oder seiner Leistung zu haben.

Wenn wir bei diesen Tests die ultimative Leistung benötigen würden, müssten wir natürlich das Beste aus dem Prozessor und dem Speicher herausholen, indem wir parallelen Code schreiben (oder eine Version unter OpenCL schreiben) ... Außerdem habe ich keine Ahnung, was passieren wird, wenn ich schreibe sehr lange Ketten von Modifikatoren. Ist es einfach, Compiler-Meldungen zu debuggen und zu lesen, wenn Bereiche in komplexen Projekten verwendet werden? Erhöht die Kompilierungszeit. Ich hoffe, dass jemals jemand über diesen Artikel schreibt.

Als ich diese Tests schrieb, wusste ich selbst nicht, was passieren würde. Jetzt weiß ich - Bereiche verdienen es definitiv, in einem realen Projekt unter den Bedingungen getestet zu werden, unter denen Sie sie verwenden möchten.

Ich ging zum Basar, um einen Samowar zu kaufen.

Nützliche Links

Range-V3 nach Hause

Dokumentation und Fallstudien v3

Code aus diesem Artikel auf Github

Listen in Haskell zum Vergleich

Danksagung

Danke fadey, dass du beim Schreiben geholfen hast !

PS

Ich hoffe, jemand kommentiert solche Kuriositäten: i) Wenn Sie das Integrationsintervall zehnmal kleiner nehmen, ist das Beispiel v2 auf meinem Xeon 10% schneller als v1 und v4 dreimal schneller als v1. ii) Der Intel-Compiler (icc 2019) erstellt in diesen Beispielen manchmal Code, der doppelt so schnell ist wie kompiliertes g ++. Ist die Vektorisierung schuld? Kann g ++ dazu gezwungen werden?