Dies ist eine alte Version des Dokuments!


Kumulative Summe

Als kumulierte Summe einer Zahlenfolge bezeichnet man die einfache Addition der einzelnen Zahlen. Dabei wird jede Zahl zur vorhergehenden Summe addiert.
Wir können die Berechnung allgemein so darstellen: si = si-1 + xi

Um diese Definition zu verdeutlichen, berechnen wir für zufällig gewählte Zahlen die kumulierten Summe:

Zahlen 3 7 12 4 -6 2 -4
Kumulierte Summe 3 10 22 26 20 22 18

Noch einmal der Rechenvorgang aufgelistet:

Index (i) Zahl (x) Summe (s) Berechnung
1 3 3 0 + 3
2 7 10 3 + 7
3 12 22 10 + 12
4 4 26 22 + 4
5 -6 20 26 - 6
6 2 22 20 + 2
7 -4 18 22 - 4

Jetzt wissen wir was die kumulative Summe ist und wie wir sie berechnen. Sehen wir uns nun ein paar Beispiele an, die wir über diese Berechnung lösen bzw. optimieren können.

Beispiel: Optimierung von Aktienkäufen

Problemstellung

Um mit Aktien einen möglichst hohen Gewinn zu erzielen, muss der Verkaufspreis einen möglichst großen Unterschied (der natürlich positiv sein muss) zum Einkaufspreis aufweisen. Diese Aufgabe ist nicht so trivial wie sie anfangs vielleicht aussieht, denn die beste Lösung ist nicht immer der niedrigste (Einkauf) und der höchste Preis (Verkauf) über den gesamten Zeitraum. Es muss dabei auch die Reihenfolge der Kurspositionen berücksichtigt werden.
Wir nehmen an, dass wir Vorhersagen zu künftigen Aktienkursen haben und den optimalen Einkaufs- und Verkaufszeitpunkt, sowie die jeweiligen Preise bestimmen wollen.

Eingabe

Als Eingabe bekommen wir den Kurs einer Aktie in einem Zeitintervall. Dieses Intervall kann beliebig groß sein und wird in diesem Beispiel vernachlässigt.

Wir haben also folgende Eingabegrößen:

Variable Wert
Anzahl an Kursen (n) 8

Danach folgen n Kurse:

Index (i) Kurs (k)
1 11
2 13
3 4
4 6
5 10
6 5
7 11
8 9

In unseren Test-Programmen werden die Daten über die Standardeingabe nach folgendem Format geliefert:

n
k_1
...
k_n

Bzw. in ausgeschriebener Form:

input.txt
8
11
13
4
6
10
5
11
9

Wie es für Aktienkurse üblich ist, können wir sie auch mit einem Diagramm darstellen:

Quadratische Lösung

Eine Möglichkeit die maximale Differenz zu finden, wäre die Berechnung der Differenzen zu allen nachfolgenden Kursen, aus denen man anschließend das Maximum auswählt.

squared_stock_prices.c
#include <stdio.h>
#include <stdlib.h>
 
 
// Maximum an eingegebenen Kurspositionen
const int MAX_POSITIONS = 1000000;
// Puffer für Eingaben
const int BUFFER_SIZE = 12;
 
 
int main()
{
  // Puffer für Eingaben
  char buffer[BUFFER_SIZE];
  // Zählvariablen für Schleifen
  int i, j;
  // Anzahl der Kurspositionen
  int n;
  // Werte der einzelnen Kursposition
  int v[MAX_POSITIONS];
  // Index des minimalen Wertes (optimaler Einkauf)
  int minIndex = -1;
  // Index des maximalen Wertes (optimaler Verkauf)
  int maxIndex = -1;
 
  // Anzahl an Kursen einlesen
  fgets(buffer, BUFFER_SIZE, stdin);
  n = atoi(buffer);
  if (n < 1 || n > MAX_POSITIONS)
  {
    printf("Anzahl der Kursposition muss zwischen 1 und %d liegen\n", MAX_POSITIONS);
    return 1;
  }
 
  // Kurse einlesen
  for( i = 0; i < n; i++ )
  {
    fgets(buffer, BUFFER_SIZE, stdin);
    v[i] = atoi(buffer);
  }
 
  // Über alle Kurse iterieren
  for( i = 0; i < n; i++ )
  {
    // Differenz zu allen darauffolgenden Kursen bilden.
    // Ist das der erste Durchlauf (minIndex == -1) oder
    // die Differenz größer als das bisherige Maximum,
    // werden die bisherigen Indizes überschrieben.
    for( j = i + 1; j < n; j++ )
    {
      if( minIndex == -1 || v[j] - v[i] > v[maxIndex] - v[minIndex] )
      {
        maxIndex = j;
        minIndex = i;
      }
    }
  }
 
  // Ergebnis ausgeben
  if( minIndex != -1 )
    printf( "Minimum bei Zeiteinheit: %d; Minimaler Wert: %d\n"
            "Maximum bei Zeiteinheit: %d; Maximaler Wert: %d\n"
            "Gewinn: %d\n", 
            minIndex + 1, v[minIndex],
            maxIndex + 1, v[maxIndex],
            v[maxIndex] - v[minIndex] );
 
  return 0;
}

Aufruf:

$ ./squared_stock_prices < input.txt

Ausgabe:

Minimum bei Zeiteinheit: 3; Minimaler Wert: 4
Maximum bei Zeiteinheit: 7; Maximaler Wert: 11
Gewinn: 7

Die optimale Strategie im angegebenen Zeitintervall ist also im 3. Zeitintervall zu kaufen und im 7. zu verkaufen. Dabei erhalten wir einen Gewinn von 7 pro Aktie.

Sehen wir uns nun die Anzahl der für die Berechnung notwendigen Schritte an. Folgende beiden Schleifen sind für die Berechnung der kumulativen Summe verantwortlich:

  // Über alle Kurse iterieren
  for( i = 0; i < n; i++ )
  {
    // Differenz zu allen darauffolgenden Kursen bilden.
    // Ist das der erste Durchlauf (minIndex == -1) oder
    // die Differenz größer als das bisherige Maximum,
    // werden die bisherigen Indizes überschrieben.
    for( j = i + 1; j < n; j++ )
    {
      if( minIndex == -1 || v[j] - v[i] > v[maxIndex] - v[minIndex] )
      {
        maxIndex = j;
        minIndex = i;
      }
    }
  }

Wir haben zwei Schleifen, die jeweils einen Zähler gegen n laufen lassen. Die äußere Schleife startet bei 0, hat also genau n Iterationen.

In der inneren Schleife starten wir bei i + 1, die Anzahl der Schritte ist somit von der äußeren Schleife abhängig. Der aufwendigste Fall mit den meisten Schritten ist der allererste, bei dem i den Wert 0 hat und die innere Schleife bei 1 startet. In diesem Fall hat die innere Schleife n - 1 Schritte. Nachdem wir nur eine grobe Schätzung wollen, lassen wir den konstanten Wert weg und runden das Ergebnis einfach auf n. Nun haben wir für jeden der n Durchläufe der äußeren Schleife auch n Durchläufe in der inneren Schleife. Es ergibt sich ein von der Anzahl an Kursen n abhängiger quadratischer Aufwand. In Kurzschreibweise bedeutet das: O(n2)

Lineare Lösung: Kumulative Summe

Die obrige Lösung funktioniert in unserem Beispiel ganz gut. Aber was ist, wenn wir nicht 8 sondern ein paar Millionen Kurspositionen haben? Versuchen wir eine bessere Lösung mit der kumulativen Summe zu finden.

Gehen wir von der ersten Kursposition aus und addieren alle Kursdifferenzen dazu, erhalten wir immer den aktuellen Kurs. Wenn wir uns jetzt das Minimum merken und immer die Differenz zum aktuellen Kurs bilden, können wir in linearem Aufwand das gesuchte Ergebnis finden.

linear_stock_prices.c
#include <stdio.h>
#include <stdlib.h>
#include <limits.h>
 
 
// Maximum an eingegebenen Kurspositionen
const int MAX_POSITIONS = 1000000;
// Puffer für Eingaben
const int BUFFER_SIZE = 12;
 
 
int main()
{
  // Puffer für Eingaben
  char buffer[BUFFER_SIZE];
  // Anzahl der Kurspositionen
  int n;
  // Index des minimalen Wertes (optimaler Einkauf)
  int minIndex = 0;
  // Index des maximalen Wertes (optimaler Verkauf)
  int maxIndex = -1;
  // Minimaler Wert (Einkaufspreis)
  int min;
  // Maximale Kursdifferenz (Verkaufspreis)
  int maxDiff = INT_MIN;
  // Aktuelle Summe aller Kursänderungen. Dadurch enthält
  // diese Variable auch immer den aktuellen Kurs.
  long sum;
  // Aktueller Aktienkurs
  int currentPrice;
  // Vorheriger Aktienkurs
  int oldPrice;
 
  // Anzahl an Kursen einlesen
  fgets(buffer, BUFFER_SIZE, stdin);
  n = atoi(buffer);
  if (n < 1 || n > MAX_POSITIONS)
  {
    printf("Anzahl der Kursposition muss zwischen 1 und %d liegen\n", MAX_POSITIONS);
    return 1;
  }
 
  // 1. Kurs einlesen und als bisherige Summe festlegen.
  fgets(buffer, BUFFER_SIZE, stdin);
  currentPrice = atoi(buffer);
  sum = min = oldPrice = currentPrice;
  // Alle Aktienkurse einlesen und Ergebnis berechnen
  int i;
  for( i = 1; i < n; i++ )
  {
    // Aktuellen Kurs einlesen
    fgets(buffer, BUFFER_SIZE, stdin);
    currentPrice = atoi(buffer);
    // Wir summieren die Kursänderungen auf.
    sum += currentPrice - oldPrice;
    // Ist der Kurs kleiner als das bisherige Minimum, wird 
    // die aktuelle Kursposition als neues Minimum gespeichert.
    if( sum < min )
    {
      min = sum;
      minIndex = i;
    }
    // Ist der Kursunterschied zwischen dem bisherigen minimalen und
    // dem aktuellen Kurs größer als das bisherige Maximum, wird
    // die aktuelle Kursposition als neues Maximum gespeichert.
    if( sum - min > maxDiff )
    {
      maxDiff = sum - min;
      maxIndex = i;
    }
    // Wir speichern den aktuellen Kurs für den nächsten Durchlauf als
    // alten Kurs.
    oldPrice = currentPrice;
  }
 
  // Ergebnis ausgeben.
  printf( "Minimum bei Zeiteinheit: %d; Minimaler Wert: %d\n"
          "Maximum bei Zeiteinheit: %d; Maximaler Wert: %d\n"
          "Gewinn: %d\n", 
          minIndex + 1, min, 
          maxIndex + 1, maxDiff + min, 
          maxDiff );
 
  return 0;
}

Aufruf:

$ ./linear_stock_prices < input.txt

Ausgabe:

Minimum bei Zeiteinheit: 3; Minimaler Wert: 4
Maximum bei Zeiteinheit: 7; Maximaler Wert: 11
Gewinn: 7

In dieser Version haben wir nicht nur den rechnerischen Aufwand von O(n2) auf O(n) verringert, wir ersparen uns sogar ein Array mit den gesamten Daten.

Performance-Vergleich

Nachfolgend wollen wir einen einfachen Performance-Vergleich zwischen den beiden Implementierungen durchführen. Es ist zu beachten, dass dies kein vollständiger und statistisch zuverlässiger Benchmark ist. Vielmehr soll der Sinn hinter dieser Art von Optimierung demonstriert werden.

Für einen sinnvollen Vergleich sind große Datenmengen erforderlich. Um die entsprechenden Eingabedaten zu generieren, wird ein Python-Skript verwendet:

generate.py
#!/usr/bin/env python3
import random
import sys
 
n = int(sys.argv[1])
lower_bound = int(sys.argv[2])
upper_bound = int(sys.argv[3])
 
print(n)
for i in range(n):
    print(random.randint(lower_bound, upper_bound))

Nun erzeugen wir 1 Million Kurspositionen zwischen 0 und 1000. Aufruf:

$ ./generate.py 1000000 0 1000 > large_input.txt

Beide Implementierungen werden auf höchster Optimierungsstufe kompiliert:

$ gcc -O3 squared_stock_prices.c -o squared_stock_prices
$ gcc -O3 linear_stock_prices.c -o linear_stock_prices

Starten wir die Berechnung für die quadratische Lösung:

$ time ./squared_stock_prices < large_input.txt
Minimum bei Zeiteinheit: 1227; Minimaler Wert: 0
Maximum bei Zeiteinheit: 5013; Maximaler Wert: 1000
Gewinn: 1000
./squared_stock_prices < large_input.txt  1557.53s user 0.18s system 99% cpu 26:03.02 total

In diesem Fall dauerte der Durchlauf über 26 Minuten. Zum Vergleich die lineare Lösung:

$ time ./linear_stock_prices < large_input.txt                
Minimum bei Zeiteinheit: 1227; Minimaler Wert: 0
Maximum bei Zeiteinheit: 5013; Maximaler Wert: 1000
Gewinn: 1000
./linear_stock_prices < large_input.txt  0.06s user 0.00s system 98% cpu 0.061 total

Die optimierte Lösung benötigte für das gleiche Ergebnis 0.06 Sekunden! Natürlich werden bei diesem Vergleich einige wesentliche Dinge vernachlässigt, wie etwa die verwendete Hardware oder Caches. Das Ergebnis ist trotzdem relativ eindeutig.

Beispiel: Größtmögliche freie Grundstücksfläche

Problemstellung

Wir haben ein Grundstück, auf das wir ein Gebäude bauen möchten. Dieses Grundstück ist jetzt aber keine komplett leere Fläche, sondern beinhaltet auch Hindernisse wie zum Beispiel Bäume. Da wir diese nur ungern fällen möchten, suchen wir eine möglichst große rechteckige Teilfläche, auf der wir unser Gebäude platzieren können, ohne ein Hindernis entfernen zu müssen. Es interessieren uns Koordinaten und Breite bzw. Höhe des Rechtecks.

Eingabe

Zuerst bekommen wir die Seitenlängen des Grundstücks als Eingabe:

Variable Wert
Breite 6
Höhe 8

Danach folgt die Eingabe der Fläche. Diese wird dabei in Form von Zeichen dargestellt. . steht für ein freies Feld, # für ein belegtes. Wir suchen nun den größten rechteckigen Block, der nur aus freien Feldern (.) besteht. Als Beispiel nehmen wir folgendes Grundstück an:

..#..#
......
......
.....#
......
#.....
......
......

Die gesamte Eingabe besteht somit aus diesen Informationen:

6 8
..#..#
......
......
.....#
......
#.....
......
......

n^6 Lösung

Ein relativ einfaches, aber rechnerisch sehr aufwendiges Verfahren ist es, alle möglichen Felder durchzuprobieren. Dabei geht man alle einzelnen Koordinaten durch und sucht für alle davon ausgehenden Breiten und Höhen das größte freie Feld. Eine Implementierung könnte so aussehen:

n6_free_area.c
#include <stdbool.h>
#include <stdio.h>
#include <stdlib.h>
 
const char OBSTACLE_CHARACTER = '#';
 
int main()
{
    // Dimensionen des Grundstücks einlesen
    int width, height;
    scanf("%d %d\n", &width, &height);
 
    // Grundstücks einlesen
    char *area = (char *)malloc(width * height * sizeof(char));
    for (int i = 0; i < height; i++)
    {
        for (int j = 0; j < width; j++)
            area[i * width + j] = getchar();
        // Newline ignorieren
        getchar();
    }
 
    // Noch keine passende Fläche gefunden, deshalb mit -1 vorbelegt
    int solution_x = -1, solution_y = -1, solution_w = -1, solution_h = -1;
    bool solution_found = false;
 
    // Durch alle Zeilen iterieren (y-Koordinate)
    for (int start_y = 0; start_y < height; start_y++)
    {
        // Durch alle Spalten iterieren (x-Koordinate)
        for (int start_x = 0; start_x < width; start_x++)
        {
            // Alle möglichen Höhen versuchen, beginnend bei der größten
            for (int end_y = start_y; end_y < height; end_y++)
            {
                // Alle möglichen Breiten versuchen, beginnend bei der größten
                for (int end_x = start_x; end_x < width; end_x++)
                {
                    // Prüfen, ob der Bereich leer ist
                    unsigned int num_window_obstacles = 0;
                    for (int window_y = start_y;
                         window_y <= end_y && num_window_obstacles == 0;
                         window_y++)
                    {
                        for (int window_x = start_x;
                             window_x <= end_x && num_window_obstacles == 0;
                             window_x++)
                        {
                            if(area[window_y * width + window_x] == OBSTACLE_CHARACTER)
                                num_window_obstacles++;
                        }
                    }
 
                    // Überprüfen, ob in dem Bereich Hindernisse gefunden wurden
                    if (num_window_obstacles == 0)
                    {
                        // Prüfen, ob das neue Feld größer ist als das bisher
                        // größte gefundene. Gegebenenfalls Koordinaten und
                        // Größe speichern
                        int candidate_w = end_x - start_x + 1;
                        int candidate_h = end_y - start_y + 1;
                        if (!solution_found || ((candidate_w * candidate_h) > (solution_w * solution_h)))
                        {
                            solution_found = true;
                            solution_x = start_x;
                            solution_y = start_y;
                            solution_w = candidate_w;
                            solution_h = candidate_h;
                        }
                    }
                    else
                    {
                        // Die Schleifen verkleinern die zu überprüfende Fläche.
                        // In der innersten Schleife kann somit keine größere Fläche
                        // mehr gefunden werden.
                        break;
                    }
                }
            }
        }
    }
 
    // Ergebnis ausgeben
    if (solution_found)
        printf("x: %d y: %d w: %d h: %d -> %d\n", solution_x,
                                                  solution_y,
                                                  solution_w,
                                                  solution_h,
                                                  solution_w * solution_h);
    else
        printf("Keine freie Flaeche gefunden");
 
    free(area);
    return 0;
}

Ausgabe:

x: 1 y: 1 w: 4 h: 7 -> 28

Das Ergebnis sieht schon mal sehr gut aus. Aber am Code fällt sofort die Verschachtelung von 6 Schleifen auf. Wir haben hier also einen Algorithmus der Ordnung O(n6), der sehr schnell anwächst.

n^4 Lösung: Kumulative Summe

So möchten wir unser Programm aber nicht einfach stehen lassen. Es liefert zwar die richtige Lösung, die Berechnung ist aber sehr aufwendig und der Code durch die sechs verschachtelten Schleifen schon recht unleserlich.

n4_free_area.c
#include <stdbool.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
 
int main()
{
    // Dimensionen der Fläche einlesen
    int width, height;
    scanf("%d %d\n", &width, &height);
 
    // Fläche einlesen
    int *sum = (int *)malloc((height + 1) * (width + 1) * sizeof(int));
    memset(sum, 0, (height + 1) * (width + 1) * sizeof(int));
 
    for (int i = 1; i <= height; i++)
    {
        for (int j = 1; j <= width; j++)
        {
            int current_value = (getchar() == '#') ? 1 : 0;
            sum[i * (width + 1) + j] = current_value + sum[(i - 1) * (width + 1) + j] + sum[i * (width + 1) + j - 1] - sum[(i - 1) * (width + 1) + j - 1];
        }
        // Newline ignorieren
        getchar();
    }
 
    // Noch keine passende Fläche gefunden, deshalb mit -1 vorbelegt
    int solution_x = -1, solution_y = -1, solution_w = -1, solution_h = -1;
    bool solution_found = false;
 
    for (int start_y = 1; start_y <= height; start_y++)
    {
        for (int start_x = 1; start_x <= width; start_x++)
        {
            for (int end_y = start_y; end_y <= height; end_y++)
            {
                for (int end_x = start_x; end_x <= width; end_x++)
                {
                    unsigned int num_window_obstacles = sum[end_y * (width + 1) + end_x] -
                                                        sum[end_y * (width + 1) + start_x - 1] -
                                                        sum[(start_y - 1) * (width + 1) + end_x] +
                                                        sum[(start_y - 1) * (width + 1) + start_x - 1];
                    if (num_window_obstacles == 0)
                    {
                        int candidate_w = end_x - start_x + 1;
                        int candidate_h = end_y - start_y + 1;
                        if (!solution_found || (candidate_w * candidate_h) > (solution_w * solution_h))
                        {
                            solution_found = true;
                            solution_x = start_x - 1;
                            solution_y = start_y - 1;
                            solution_w = candidate_w;
                            solution_h = candidate_h;
                        }
                    }
                    else
                    {
                        // Die Schleifen verkleinern die zu überprüfende Fläche.
                        // In der innersten Schleife kann somit keine größere Fläche
                        // mehr gefunden werden.
                        break;
                    }
                }
            }
        }
    }
 
    // Ergebnis ausgeben
    if (solution_found)
        printf("x: %d y: %d w: %d h: %d -> %d\n", solution_x,
                                                  solution_y,
                                                  solution_w,
                                                  solution_h,
                                                  solution_w * solution_h);
    else
        printf("Keine freie Flaeche gefunden");
 
    return 0;
}



FIXME skizzen, erklärung

Performance-Vergleich

Wie bereits zuvor, möchten wir die Performance der beiden Implementierungen vergleichen und erstellen ein Programm zur zufälligen Generierung von Eingabedaten.

2d_cumulative_sum_generator.c
#include <stdbool.h>
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
 
// Bestimmt wie wahrscheinlich ein Hindernis auftritt.
// Je größer die Zahl, desto unwahrscheinlicher treten
// Hindernisse auf.
const unsigned short OBSTACLE_FACTOR = 20;
const char OBSTACLE_SYMBOL = '#';
const char FREE_SYMBOL = '.';
 
int main(int argc, char *argv[])
{
    if (argc != 3)
    {
        printf("Usage: %s width height\n", argv[0]);
        return 1;
    }
 
    // Der Zufallsgenerator wird mit der aktuellen Zeit initialisiert.
    srand(time(NULL));
 
    // Breite und Höhe des zu generierenden Feldes werden als
    // Kommandozeilenparameter übergeben.
    int width = atoi(argv[1]),
        height = atoi(argv[2]);
 
    // Dimensionen des Grundstücks ausgeben.
    printf("%d %d\n", width, height);
    // In den Schleifen folgt die zufällige Generirung
    // und Ausgabe des Grundstücks.
    for (int i = 0; i < height; i++)
    {
        for (int j = 0; j < width; j++)
        {
            // Eine Zufallszahl wird generiert. Ist sie durch
            // den konstanten Faktor teilbar, befindet sich
            // an dieser Stelle ein Hindernis.
            bool is_obstacle = rand() % OBSTACLE_FACTOR == 0;
            if (is_obstacle)
                printf("%c", OBSTACLE_SYMBOL);
            else
                printf("%c", FREE_SYMBOL);
        }
        printf("\n");
    }
 
    return 0;
}

Das Programm bekommt die gewünschte Flächengröße als Kommandozeilenparameter und schreibt die generierten Daten auf die Standardausgabe. Für reproduzierbare Testläufe empfiehlt es sich deshalb, die Ergebnisse in einer Datei zwischenzuspeichern:

$ ./2d_cumulative_sum_generator 6 8 > 6x8.txt

Generieren wir nun ein paar größere Datensätze mit quadratischen Dimensionen, um die Performance besser einschätzen zu können.

$ for s in `seq 500 500 3000`; do ./2d_cumulative_sum_generator "$s" "$s" > "${s}x${s}.txt"; done
Datei Dauer n6 Ausführung Dauer n4 Ausführung
500x500.txt 1,318 0,186
1000x1000.txt 9,814 1,359
1500x1500.txt 31,588 5,906
2000x2000.txt 1:16,49 16,988
2500x2500.txt 2:31,85 32,424
3000x3000.txt 4:17,61 1:01,30

Es handelt sich bei diesen Tests um einzelne Ausführungen mit vorangestellter Zeitmessung über time. Abermals geht es hier nicht darum eine möglichst genaue Messung zu schaffen, sondern die praktische Größenordnung des Unterschieds einschätzen zu können.

Die Performance der n6 Implementierung hängt auch deutlich von der Häufigkeit der Hindernisse ab. Die strengen Bedinungen in den innersten beiden Schleifen stoppen beim ersten gefundenen Hindernis. Weiters wird auch die Vergrößerung der zweiten Dimension abgebrochen, wenn ein Hindernis gefunden wurde. Treten Hindernisse häufig auf, ergeben sich dadurch weniger Schleifendurchläufe. In Extremfällen kann sich die Performance sogar an jene der n4 Implementierung annähern.

Siehe auch