Ein Metaheuristikframework mit Hilfe von abstrakten Klassen

Blaupause für ein simples Metaheuristik-Framework.

Ausgangslage

In diesem Artikel entwerfe ich ein simples metaheuristisches Framework und mache mir dabei abstrakte Klassen zunutze. Dadurch können verschiedene Probleme mit geringem Adaptierungsaufwand implementiert und gelöst werden.

Vorwissen: Abstrakte Klassen

Abstrakte Klassen definieren einen grundlegenden Klassenbauplan, bestehend auf Funktionen und (optional) Klassenvariablen. Hierbei unterscheiden sich abstrakte Klassen von herkömmlichen Klassen dadurch, dass sie

1. nicht selber instanziert werden können
2. ”leere” Funktionen implementieren kann, die dann für den Anwendungsfall mit Leben gefüllt werden sollen
3. Funktionen/Variablen enthalten können, die für alle Klassen, die die abstrakte Klasse implementieren, gleich sein.

Ausführlicher auf Wikipedia

Vorwissen: Simulated Annealing

Simulated Annealing ist ein zufallsbasiertes metaheuristisches Verfahren, das iterativ Lösungen entwift und verändert. Hierbei wird versucht durch sukzessiv abnehmende Toleranz von verschlechternden Lösungsveränderungen in einen Bereich des Lösungsraums vorzudringen, in dem (hoffentlich) gute bzw. optimale Lösungen zu finden sind.

Ausführlicher auf Wikipedia

Herangehensweise

Egal welches Problem mit Simulated Annealing (oder anderen metaheuristischen Frameworks) gelöst werden soll, der Ansatz ist immer sehr ähnlich:

Eine Startlösung wird kleinschrittig verändert, in der Hoffnung sich so einer optimalen Lösung zu nähern.

Beispiele:

• Bei der Tourenplanung werden z.B. Knoten in der Tour getauscht
• Bei der Standortplanung werden Standorte geöffnet oder geschlossen
• Bei der Produktionsplanung werden Rüstzeitpunkte festgelegt

Wir müssen also sicherstellen, dass zwei zentrale Konzepte implementiert werden:

• Ein Datenkonstrukt, das eine Lösung zum betrachteten Problem abbildet
• Eine Methode, das die Lösung verändern kann

Die konkreten Ausprägungen der beiden Konzepte sind stark problemspezifisch. Die Lösung zu einem Tourenplanungproblem ist i.d.R. eine Sequenz von besuchten Orten, während die Lösung zu einem Standortplanungsproblem meist eine Kombination aus eröffneten Standorten und Materialflüssen darstellt.

Abstrakte Lösungs- und Manipulatorklasse

Nachfolgend werden zunächst die beiden notwendigen abstrakten Klassen entworfen.

Abstrakte Lösungsklasse

Die Klasse wird so knapp wie möglich geschrieben. Für das Verfahren notwendig, ist lediglich der Zielfunkionswert objective. Ich neige darüber hinaus dazu, Lösungen immer erst zu kopieren und anschließend zu verändern.

Deswegen implemetiere ich eine copy() Methode, die mir eine DeepCopy (Link) der Lösungsinstanz zurück gibt. Wenn ich diese Kopie nachfolgend verändere, bleibt die Ursprungsinstanz unberührt.

from abc import ABC, abstractmethod

class AbstSolution(ABC):

    """
        Abstrakte Lösungsklasse. 
        Alle Implementierungen enthalten:
        - Zielfunktionswert der Lösung
        - Funktion, die eine DeepCopy der aktuellen Klasseninstanz erzeugt
    """

    def __init__(self):
        # Nur der Zielfunktionswert ist vorgegeben
        self.objective = 0
    
    @abstractmethod
    def copy(self):
        # Soll eine DeepCopy erzeugen
        pass

Abstrakte Lösungsveränderungsklasse

Die Lösungsveränderungsklasse bleibt ebenfalls simpel: Hier wird nur die Funktion genNextSol angelegt, mit der der nächste Lösungskandidat generiert werden soll.

from abc import ABC, abstractmethod
class AbstSolutionMnpltr(ABC):

    @abstractmethod
    def genNextSol(self, incsol : AbstSolution) -> AbstSolution:
        """
        Manipuliert eine Lösung + gibt sie zurück.
        """
        pass

Simulated Annealing auf Basis der abstrakten Klassenkonstrukte

Nachfolgend wird das Simulated Annealing Verfahren implementiert. Wir stellen fest: an keiner Stelle basiert die Implementierung auf einer konkreten Ausgestaltung der abstrakten Klassen. Es wird lediglich vorausgesetzt, dass irgendwie Lösungen verändert werden (nextSolution = solManipulator.genNextSol(incSol)).

In __isBetter wird dann nur anhand des objectives geprüft, ob Lösungen besser oder schlechter sind.

Es liegt also am Anwender, die Klassen mit problemspezifischen Verfahren zu füllen.

class SimulatedAnnealing:

    def __init__(self, minimize : bool):
        # Gibt an, ob minimiert oder maximiert werden soll.
        self.minimize = minimize
        self.bestSolution : AbstSolution
    
    def __isBetter(self,solA: AbstSolution, solB : AbstSolution):
        """ Abh. davon ob minimiert oder maximiert wird: prüft
        anhand des Objectives ob eine Lösung A besser ist als 
        eine Lösung B"""
        if self.minimize:
            return solA.objective < solB.objective
        
        return solA.objective > solB.objective
    
    def __calcExp(self, nextSol : AbstSolution, incSol : AbstSolution, t : float):
        # e^(-delta/T): nahe 1, bei kleinem Delta bzw. großem T, nahe 0 sonst.
        delta = abs(nextSol.objective-incSol.objective)
        return exp((-(delta))/t) 

    def solve(self, startSol : AbstSolution, solManipulator: AbstSolutionMnpltr, startTemp : int, stopTemp : int, alpha: float):

        # Initialisiert beste Lösung
        self.bestSolution  = startSol.copy()
        incSol = startSol

        t = startTemp

        while(t > stopTemp):

            # Erzeugt veränderte Lösung
            nextSolution = solManipulator.genNextSol(incSol)

            if self.__isBetter(nextSolution, incSol):  
                # Falls besser als aktuelle Lösung:
                # aktualisierte aktuelle Lösung
                incSol = nextSolution
                if self.__isBetter(nextSolution,self.bestSolution):    
                    # Falls besser als beste bekannte Lösung:
                    # aktualisiere beste Lösung
                    self.bestSolution = nextSolution.copy()
                    print(self.bestSolution.objective)
            else:
                # Falls nicht besser: 
                # akzeptiere nur mit gewisser Wahrscheinlichkeit
                if random.random() <= self.__calcExp(nextSolution,incSol, t): 
                    incSol = nextSolution

            t = alpha * t  

Konkrete Implementierung: Tourenplanung

Die beiden abstrakten Klassen werden jetzt für ein Tourenplanungsproblem implementiert. Wir nehmen an, dass es eine Anzahl von Orten gibt, von denen jeder im Rahmen einer einzelnen Rundtour exakt einmal besucht werden soll und dass diese Tour an einem Ort “0” startet und endet.

Hierfür haben wir eine Distanzmatrix mit Distanzen zwischen jedem Paar von Orten gegeben. Der Ort an Index 0 stellt den Start- bzw. Zielort dar.

Tour

Die Lösungsklasse implementieren wir, indem wir die Sequenz und den Zielfunktionwert speichern.

class Tour(AbstSolution):
    def __init__(self):
        super().__init__()
        self.objective = 0 # Länge der Tour
        self.sequence = [] # Sequenz besuchter Orte

    def copy(self) -> AbstSolution:
        # Kopiert die aktuelle Lösung
        s = Tour()
        s.objective = self.objective
        s.sequence = self.sequence[:]
        return s 

Tourveränderung

Nachfolgend wird eine Klasse entworfen, die den AbstSolutionMnpltr implementiert. In dieser werden zufällig zwei Indizes in der Tour ausgewählt und die entsprechenden Orte werden getauscht.

Hierbei ist genNextSol die einzige Funktion, die wir durch die Erweiterung der abstrakten Klasse implementieren müssen. Die anderen Funktionen sind spezifisch für den Anwendungsfall entworfen worden.

Anstatt die Länge der Tour komplett neu zu berechnen, ermitteln wir die Veränderung der Länge (__swapDelta) und addieren diese auf die Tourlänge. Verkürzt sich die Tour, ist der Rückgabewert von __swapDelta negativ.

from abc import ABC
from random import sample
class SwapInsert(AbstSolutionMnpltr):

    def __init__(self, dist) -> None:
        super().__init__()
        self.dist = dist
    
    def __swapDelta(self, a : int, b : int, sequence):
        """
        Berechnet das Kostendelta, wenn index a 
        und index b in Sequenz getauscht werden.
        """
        # Macht es einfacher: at ist kleinerer Index
        at = min(a,b)
        bt = max(a,b)

        nodeA = sequence[at]
        aPred = self.__getPred(at,sequence)

        nodeB = sequence[bt]
        bSucc = self.__getSucc(bt,sequence)

        if at != bt-1:
            # Fall 1: a ist nicht direkter Vorgänger von b
            aSucc = self.__getSucc(at,sequence)
            bPred = self.__getPred(bt,sequence)
            delta = -(self.dist[aPred][nodeA] + self.dist[nodeA][aSucc])
            delta -= (self.dist[bPred][nodeB] + self.dist[nodeB][bSucc])
            delta += (self.dist[aPred][nodeB] + self.dist[nodeB][aSucc])
            delta += (self.dist[bPred][nodeA] + self.dist[nodeA][bSucc])
            return delta

        # Fall 2: a ist direkter Vorgänger von b
        delta = -(self.dist[aPred][nodeA] + self.dist[nodeB][bSucc])
        delta -= self.dist[nodeA][nodeB]
        delta += (self.dist[aPred][nodeB] + self.dist[nodeA][bSucc])
        delta += self.dist[nodeB][nodeA]
        return delta

    def __getPred(self, idx,sequence):
        """ Gibt den Vorgänger vom index idx in der Sequenz zurück """
        return 0 if idx == 0 else sequence[idx-1]

    def __getSucc(self,idx,sequence):
        """ Gibt den Nachfolger vom index idx in der Sequenz zurück """
        return 0 if idx == len(sequence)-1 else sequence[idx+1]

    def genNextSol(self, incsol: AbstSolution):
        """Erzeugt eine neue Lösung indem 
        zwei Orte in der Sequenz getauscht werden."""

        # Wähle zwei zufällige Knoten in der Tour
        a,b = random.sample(range(len(incsol.sequence)), 2)

        # Berechne Kostenveränderung
        costDelta = self.__swapDelta(min(a,b),max(a,b),incsol.sequence)

        # Gibt neue Lösung zurück
        newSol = incsol.copy()
        newSol.objective += costDelta
        newSol.sequence[a], newSol.sequence[b] = newSol.sequence[b], newSol.sequence[a]
        return newSol

Testlauf

Zu guter letzt testen wir das Framework anhand einer zufälligen Instanz mit 100 Orten und symmetrischen Distanzen. Hierbei starten wir mit der Tour 1,2,3,4....,NR_NODES.

NR_NODES = 100
dist = [] # Distanzmatrix
for n in range(NR_NODES+1):
    dist.append([0]*(NR_NODES+1))
    for i in range(0,n):
        dist[i][n] = random.randint(100,200)
        dist[n][i] = dist[i][n]

# Initialisiere die Klassen
startSol = Tour()
solGenerator = SwapInsert(dist=dist)

# Startlösung [1,2,3,4,...,NR_NODES]
startSol.sequence = list(range(1,NR_NODES+1))
startSol.objective = dist[0][1] + dist[NR_NODES][0]
for i,n in enumerate(startSol.sequence):
    if i==0:
        startSol.objective += dist[0][n]
    elif i==len(startSol.sequence)-1:
        startSol.objective += dist[n][0]
    else:
        startSol.objective += dist[startSol.sequence[i-1]][n]

simAnneal = SimulatedAnnealing(minimize=True)
simAnneal.solve(startSol,solGenerator,startSol.objective*1.5,1,0.99)

Try yourself

Der oben dargestellte Code liegt auf GitHub.

Ausprobiert werden kann er auf Google Colab