Allokationsmodul für Voith-Schneider und Azimutpropeller

Einleitung

Das Allokationsmodul entstand im Rahmen eines Promotionsvorhaben, weshalb sich diese Dokumentation in einigen Teilen auf die Dissertation1 bezieht. Das Allokationsmodul ist darauf ausgelegt, eine Steuerung für die am Schiff vorhandenen Antriebe zu finden, sodass eine von einem Regler oder einem Masterjoystick geforderte Kraft so gut wie möglich und energieeffizient erzeugt wird. Des Weiteren können nach Belieben Interaktionen, wie Propeller-Strömung-Interaktion, Propeller-Propeller-Interaktionen und Propeller-Schiff-Interaktion bei der Berechnung der Steuerung einbezogen werden.

Theoretische Grundlagen

Mathematisch lässt sich das Allokationsproblem als Optimalsteuerungsproblem mit multikriterieller Zielfunktion beschreiben (vgl. Dissertation1, Kapitel 3):

Minimiere

\[\begin{equation} \left ( \int_0^T \left \Vert \hat{\tau} - \tau_{thr} (z(t)) \right \Vert^2 dt,\ \left \Vert P(z(t))-P_{\min} \right \Vert^2 dt \right ), \end{equation} \]

sodass die Nebenbedingungen für \(t\in [0, T]\)

\[ \begin{align} \dot z(t)&= \tilde {z} \omega(t) \\ z_{\min} & \leq z(t) \leq z_{\max} \\ M_{\min} & \leq M(z(t)) \leq M_{\max} \\ P_{\min} & \leq P(z(t)) \leq P_{\max} \\ \Vert \omega(t) \Vert_{\infty} & \leq 1 \\ z(0) & = z_0 \end{align} \]

erfüllt sind. Nach einer Diskretisierung des Problems erhält man ein nichtlineares Optimierungsproblem mit Nebenbedingungen.

Methoden der multikriteriellen Optimierung

Allgemein lässt sich ein multikriterielles bzw. ein bikriterielles Optimimierungsproblem wie folgt formulieren:

\[ \begin{equation} f(x):= (f_1(x),f_2(x)) \rightarrow \min \end{equation} \]

Es gibt verschiedene Methoden, diese Probleme zulösen. Die beiden im Modul verwendeten Methoden werden nun kurz vorgestellt. Einzelheiten lassen sich in der Dissertation1 in Kapitel 3 nachlesen.

Lexikographische Optimierung

Die beiden Zielfunktionen werden hier nacheinander optimiert, d.h. :

Erste Optimierung:

\[\begin{equation} \tilde{x}= \min\ f_1(x) \end{equation} \]

Sei nun \(\mathcal{X_1}\) die Menge der Lösungen der ersten Optimierung.

Zweite Optimierung:

\[\begin{equation} \hat{x}= min\ f_2(x),\ \hat{x}\in \mathcal{X}_1 \end{equation} \]

Elastic-Constraint Methode (Modifikation)

Bei der Elastic-Constraint Methode wird eine Zielfunktion minimiert und die andere als “Strafterm” hinzugefügt.

\[\begin{equation} \int_0^T\left\Vert P(z(t))-P_{\min}\right\Vert^2 dt +\mu\int_0^T(\left \Vert\hat{\tau}-\tau_{thr}(z(t)) \right \Vert^2-\varepsilon)_+ dt \longrightarrow min \end{equation} \]

Software

Zur Lösung des nichtlinearen Optimierungsproblems wird Ipopt, ein Innere-Punkte-Löser, von A. Wächter und L. T. Biegler verwendet. Dieses Programm kann unter https://projects.coin-or.org/Ipopt heruntergeladen werden. Dort befindet sich auch eine Dokumentation zur Installation des Programms. Das Allokationsmodul wurde mit der Version 11.3.5 getestet.

Inhalt des Allokationsmoduls

Propellertypen

Methoden der multikriteriellen Optimierung

Diskretisierungen

Interaktionen

Struktur des Moduls

Konfigurationsdatei

Die Konfigurationsdatei enthält alle für die Allokation wichtigen Informationen. Diese Informationen ändern sich während der Berechnungen im Gegensatz zur Propellersteuerung und der geforderten Kraft nicht. Diese Datei befindet sich im Ordner “Konfigurationen” und wird als Textdatei, Datei.txt, gespeichert.

Allgemeine Informationen

Eingangssignal: 1
Ausgangssignal: 1
Wiederholungen: 20
DP: 0

Das Eingangssignal gibt an, in welcher Häufigkeit eine neue Kraftforderung an das Modul gesendet wird. In der theoretischen Formulierung entspricht das Eingangssignal dem Endzeitpunkt \(T\) im Optimalsteuerungsproblem. Dieses gibt an, wie oft ein neues Steuersignal an die Propeller gesendet wird. In der mathematischen Formulierung entspricht dies der Diskretsierung der Differenzialgleichung. Sowohl das Eingangssignal als auch das Ausgangssignal werden in der Einheit Sekunden angegeben. Für den Fall, dass Parameterstudien durchgeführt werden sollen, bei denen über mehrere Zeitintervalle die selbe Kraft gefordert wird, kann bei Wiederholungen eingestellt werden, wie oft die Allokation mit dieser Kraft durchgeführt werden soll. Hiermit kann beispielsweise untersucht werden, wie sich die Antriebe bei konstanten Kräften verhalten. Zuletzt wird bei DP eingestellt, ob die geforderte Kraft skaliert wird oder nicht. Skalierungen sind vor allem dann notwendig, wenn die geforderte Kraft nicht durch die Antriebe erzeugt werden kann und eine zeitoptimale Steuerung verwendet werden soll.

Zusammenfassend gilt:

Anzahl der Propeller

In diesem Block wird angegeben, wie viele Propeller jedes Anriebstyps verwendet werden sollen. Die Anzahl der einzelnen Propeller ist dabei eine natürliche Zahl.

Anzahl_VSP: 2
Anzahl_Azimuth: 1
Anzahl_Bugstrahler: 1

Interaktionen

Eine Interaktion wird mit 1 aktiv und mit 0 nicht beachtet. In diesem Fall wird nur die Propeller-Propeller-Interaktion aktiv.

Interaktion_Anstroemung: 0
Interaktion_Propeller_Propeller: 1
Interaktion_Propeller_Schiff: 0

Zielfunktion der multikriteriellen Optimierung und Diskretisierung

Bei der Zielfunktion kann zwischen der lexikographischen Optimierung, der zeitoptimalen Optimierung und der Elastic-Constraint Methode gewählt werden

Zielfunktion: Lexiko/Lexiko_zeit/Elastic

Die Gewichte

Gewicht_x: 1.e0
Gewicht_y: 1.e0
Gewicht_Mz: 1.e0

gewichten die drei Kraftkomponenten der geforderten Kraft in der Zielfunktion.

Anmerkung: Bei der lexikographischen Optimierung sollten die Gewichte zwischen 0 und 10 gewählt werden. Bei der Elastic-Constraint Methode müssen sie aufgrund der Leistung, die ebenfalls Einfluss auf die Zielfunktion hat, im Bereich von 1.e10 gewählt werden.

Die Toleranz hat je nach ausgewählter Zielfunktion eine unterschiedliche Bedeutung.

\[(\tau_{thr_1}(z(t))-\tau_{thr}(z(t)))^2\leq tol \]
Toleranz: 1.e-6

Für das Diskretisierungsverfahren der Diffenzialgleichung kann wahlweise zwischen dem explizitem Eulerverfahren und dem klassischen Runge-Kutta-Verfahren (vierter Ordnung) gewählt werden.

Diskretisierung: Euler/Runge

Propellerbeschreibung

Um die verschiedenen Propeller unterscheiden zu können, werden die Durchmesser, minimalen Drehzahlen usw. der Antriebe mit dem Propellernamen kombiniert, z.B.

Name: Azi1
Durchmesser_Azi1: 3

Jeder Name eines VSP beginnt mit VSP, eines Azimutpropellers mit Azi und eines Bugstrahlers mit Bug.

Voith-Schneider-Propeller

Azimutpropeller

Bugstrahlruder

Bei einer potenziellen Anströmung wird der Winkel \(\beta\) zwischen Anströmung und Propellerstrahl betrachtet.

Anmerkung: Beginnt eine Zeile mit einer Raute #, so wird sie als Kommentar gewertet.

Eingabedaten

Die Eingabedaten haben den selben Dateinamen wie die Konfigurationsdatei und enden mit .dat, Datei.dat. Während des DP-Betriebes besteht diese Datei aus einer Zeile, die mit jeder neuen Kraftanforderung aktualisiert wird. Für Fall- und Parameterstudien ist es möglich, mehrere Zeilen in diese Datei zu schreiben und damit unterschiedliche Kraftanforderungen direkt zu untersuchen. Die Einträge in der Zeile werden durch Strichpunkte “;” getrennt. Beendet wird eine Zeile ebenfalls durch einen Strichpunkt.

Die ersten drei Einträge sind die geforderten Kräfte in x- und y-Richtung sowie das Moment um die z-Achse; jeweils in Newton bzw. Newtonmeter. Dann folgt die aktuelle Steuerung der Propeller: (x-Steuerung aller VSP \(\in[-0.8, 0.8]\), y-Steigung aller VSP \(\in[-0.8, 0.8]\), Drehzahl aller VSP in [rps], Azimutwinkel aller Azimutpropeller in [rad], Drehzahl aller Azimutpropeller in [rps], Steuerung aller Bugstrahler \(\in [-1, 1]\). Im Anschluss folgt für jeden einzelnen Propeller die Anströmung in x- und y Richtung, angegeben in [m/s]. Bei einer aktiven Interaktion zwischen Propeller und Schiff folgt zuletzt für jeden einzelnen Propeller, außer für die Busgtrahler, der prozentuale Restschub bei Schub in diskrete Richtungen. Dabei ist darauf zu achten, dass die diskreten Schubrichtungen, in denen die Schubreduktion gemessen wird, äquvidistante Abstände haben.

Folgendes Beispiel zeigt eine Eingabe mit einem VSP, einem Azimutpropeller und einem Bugstrahler mit einer aktiven Propeller-Schiff-Interaktion:

0; 0; -8000000; 0.2; 0.3; 0.6; 0; 1; 0.5; 0; 0; -1; -1; 2; 0; 80; 85; 90; 85; 80; 85; 90; 85;

Es wird also \(\hat{\tau}=(0,0,-8000)^T\) kN gefordert. Die Startsteuerung der Antriebe ist für

Der VSP wird nicht angeströmt, der Azimutpropeller mit (-1, -1) m/s und der Bugstrahler mit (2, 0) m/s. Insgesamt folgen nun acht weitere Einträge für die Schubreduktion bei der Propeller-Schiff-Interaktion. Da es insgesamt zwei Hauptantriebe (VSP, Azimutpropeller) gibt, werden jedem Antrieb vier Einträge zugeteilt. Aufgrund der Wahl von äquvidistanten Schubrichtungen gilt für beide Antriebe:

Ergebnisse

Die Ergebnisse werden in zwei Dateien im Ordner “Ergebnisse” gespeichert. Die Datei Ergebnis_Datei_Auswertung.dat enthält Informationen wie Propellerschub und Leistung aller Antriebe.

Die Lösung des Optimalsteuerungsproblems, also die Änderung der Propellersteuerung, findet sich in der Datei Ergebnis_Datei_Steuerung.dat, wobei sich die Ergebniszeile wie folgt aufbaut:

Zusätzliche Daten

Die Modelle für die Schub- und Leistungsberechnung im Allokationsmodul basieren auf CFD-Daten. Für einen flexiblen Datentausch im Modul stehen m-Files im Ordner Daten_Konfigurationen_Ipopt zur Verfügung, welche neue CFD-Daten in die für das Modul notwendige Form bringen.

Dabei ist der Speicher so gewählt, dass bei einer Aktualisierung die Daten sofort im Allokationsmodul geändert werden. Eine Kompilierung ist nach einer Datenaktualisierung nicht notwendig.

Voith-Schneider Propeller

Allgemein:

Anströmung:

Propeller-Propeller Interaktion:

Azimut: Der Azimutpropeller unterscheidet sich vom VSP bei den Interaktionen in diesem Modul nur bei der Interaktion Anströmung-Propeller. Die Daten hier sind Fourierkoeffizienten basisierend auf Messungen.

Diese Daten bestehen aus 12 Spalten mit je 21 Zeilen. Die ersten beiden Spalten sind die Koeffizienten für den Fortschrittsgrad J=0, die Spalten drei und vier enthalten die Koeffizienten für den Fortschrittsgrad 0.2. Dieses Schema wird analog für die Fortschrittsgrade J =0.4, 0.6, 0.8, 1, fortgesetzt.

Anwendung des Allokationsmoduls

Datei

Falls für Fallstudien mehrere Konfigurationen verwendet werden, können weitere Dateien in jeweils einer neuen Zeile hinzugefügt werden Diese Konfigurationen werden dann von dem Modul nacheinander ausgewertet.

Datei1
Datei2
Datei3
./Allokation_Ipopt

gestartet werden.

Implementierung

Allokation_main.cpp

Einlesen der Dateinamen der Konfigurationsdatei und der Eingabedaten

  ifstream Infos("Informationen.txt",ios::in);

Einlesen dieser Konfigurationsdateien und Eingabedaten

  Propeller_Informationen=Einlesen_Informationen::Einlesen_Propellerinformationen(Datei);

Überprüfung der Konfigurationsdaten auf Korrektheit mit Fehlercheck Konfig(...) , z.B. \(n_{\min} \leq n \leq n_{\max}\). Falls die Konfigurationsdaten nicht korrekt sind, bricht das Programm ab.

Die für die Interaktion nötigen Daten werden nur dann eingelesen, falls Interaktionen beachtet werden.

  if(Propeller_Informationen.Interaktion_bool){
    Daten=Interaktionen::Einlesen_Propellerdaten();
  }

Initialisierung der Eingabewerte: geforderte Kraft (unskaliert/skaliert)

      if(Propeller_Informationen.DP==0){    //Keine Aenderung der Eingabekraft
        tau_hat[0]=Propeller_Informationen.tau1_hat[k];    
        tau_hat[1]=Propeller_Informationen.tau2_hat[k];
        tau_hat[2]=Propeller_Informationen.tau3_hat[k];
      }
      else if(Propeller_Informationen.DP==1){    //Reduktion der gesamten Kraft
        tau_hat[0]=Propeller_Informationen.tau1_hat[k];    
        tau_hat[1]=Propeller_Informationen.tau2_hat[k];
        tau_hat[2]=Propeller_Informationen.tau3_hat[k];

Initialisierung der Eingabewerte - Startsteuerung: Verwendung der Startsteuerung aus der Eingabedatei, bzw. der resultierenden Steuerung aus dem vorherigen Durchlauf.

//VSP
      for (int l=0; l<Propeller_Informationen.m_VSP;l++){
        x_pitch_init[l]=Startkonfiguration[l];
        y_pitch_init[l]=Startkonfiguration[Propeller_Informationen.m_VSP+l];
        n_VSP_init[l]=Startkonfiguration[2*Propeller_Informationen.m_VSP+l];
      }

Überprüfung auf Zulässigkeit des Startwertes. Die Verwendung von unzulässigen Startwerten führt zu einem Abbruch des Programms.

    if((*Test_Startwert)[n1*floor(Propeller_Informationen.T/h1)]<0|(*Test_Startwert)[n1*floor(Propeller_Informationen.T/h1)]>1){
       cout <<"Achtung: Unzulaessiger Startwert!" <<endl;

Lösen des Optimierungsproblems

  status = app->OptimizeTNLP(mynlp);
  if (status == Solve_Succeeded) {
    std::cout << std::endl << std::endl << "*** The problem solved!" << std::endl;
  }
  else {
    std::cout << std::endl << std::endl << "*** The problem FAILED!" << std::endl;
  }

Zuletzt werden die Ergebnisse der Allokation, sowie die Kräfte und die Leistung der einzelnen Antriebe gespeichert.

Allokation_nlp

In Allokation_nlp.cpp wird das Optimierungsproblem aufgestellt.

Methode : get_nlp_info

bool Allokation_nlp::get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
                             Index& nnz_h_lag, IndexStyleEnum& index_style)

Anzahl der Variablen des Optimierungsproblems:

n=(*DMO).Anzahl_Variablen();

Die Anzahl der Optimierungsvariablen berechnet sich in Abhängigkeit von der Diskretisierungsgenauigkeit z=T/h und der Anzahl der einzelnen Propellertypen:

int n=z*(3*m_VSP+2*m_Azi+1*m_Bug);

Eine Ausnahme bildet die zeitoptimale Optimierung in Euler_lexiko_zeitoptimal.cpp

int n=z*(3*m_VSP+2*m_Azi+1*m_Bug)+1;

wobei \(z=5\) hier fix ist. Die zusätzliche Variable ist dabei der variable Endzeitpunkt.

Anzahl der Nebenbedingungen des Optimierungsproblems:

m = (*DMO).Anzahl_Nebenbedingungen();

Die Anzahl der Nebenbedingungen ist neben der Diskretisierungsgenauigkeit und Anzahl der Antriebe auch abhängig von der jeweiligen Methode zur Lösung des multikriteriellen Optimierungsproblems.

int m=z*(4*m_VSP+3*m_Azi+2*m_Bug);
int m=z*(4*m_VSP+3*m_Azi+2*m_Bug+3);

In diesem Fall muss neben den ursprünglichen Restriktionen, wie Leistungsbeschränkungen usw. auch der Schub zu jedem diskreten Zeitpunkt aus dem ersten Optimierungsproblem beachtet werden.

int m=z*(4*m_VSP+3*m_Azi+2*m_Bug)+3;

Die drei zusätzlichen Nebenbedingungen ergeben sich aus der Restriktion, dass die Propeller zum Endzeitpunkt die geforderte Kraft erzeugen müssen.

Anzahl der Nicht-Null-Einträge in der Jacobi-Matrix:

 nnz_jac_g = (*DMO).Anzahl_Nicht_Null_Jacobi();

Wie auch die Anzahl der Variablen und der Nebenbedingungen, hängt auch die Anzahl der Nicht-Null-Einträge der Jacobimatrix von der Art und der Anzahl der Antriebe, der Diskretisierungsgenauigkeit und der Methode der multikriteriellen Optimierung ab. Zudem ergeben sich Unterschiede bei Beachtung von Interaktionen und Nichtbeachtung derer.

Methode: get_bounds_info

bool Allokation_nlp::get_bounds_info(Index n, Number* x_l, Number* x_u,
                                Index m, Number* g_l, Number* g_u)

Beschränkung der Optimierungsvariablen:

  vector<double>* Variable_low=(*DMO).Grenze_low_Variablen();
  vector<double>* Variable_up=(*DMO).Grenze_up_Variablen();

Die Optimierungsvariablen sind beschränkt auf dem Intervall \([-1, 1]\)

for (int i=0; i<n;i++){
    (*Grenze_low_Var)[i]=-1;
  }
  for (int i=0; i<n;i++){
    (*Grenze_up_Var)[i]=1;
  }

Die Ausnahme hierzu bildet wieder die zeitoptimale Optimierung mit der Beschränkung der Variable, die den Endzeitpunkt beschreibt:

  (*Grenze_low_Var)[n-1]=0;
 
  (*Grenze_up_Var)[n-1]=1000;

Beschränkung der Nebenbedingungen:

Aufgrund einer Transformation der Nebenbedingungen \(g\) von \(g_l\leq g\leq g_u\) auf \(0\leq \frac{(g-g_l)}{(g_u-g_l)}\leq 1\) sind die Grenzen fix bei 0 und 1.

  for (Index i=0; i<m;i++){
    g_l[i]=0;
    g_u[i]=1;
  }

Methode: get_starting_point

bool Allokation_nlp::get_starting_point(Index n, bool init_x, Number* x,
                                   bool init_z, Number* z_L, Number* z_U,
                                   Index m, bool init_lambda,
                                   Number* lambda)
  vector<double>* Initial=(*DMO).Initial_Variablen();

Der Startwert für die Optimierung hängt von dem jeweiligen Optimierungsproblems ab und wird so gewählt, dass er bis auf den Sonderfall der zeitoptimalen Optimierung sicher zulässig ist.

for (int i=0; i<n;i++){
    (*Initial_Var)[i]=0;
  }

Aufgrund der zulässigen Startsteuerung ist dieser Startwert immer zulässig.

Zwischenergebnisse/Steuerung.dat gespeichert und als Initialwert für die zweite Optimierung verwendet.

  for (int i=0; i<n-1;i++){
    (*Initial_Var)[i]=0;
  }

  (*Initial_Var)[n-1]=1;  

Methode: eval_f

bool Allokation_nlp::eval_f(Index n, const Number* x, bool new_x, Number& obj_value)
  obj_value=(*DMO).Zielfunktion(x);

Die Zielfunktion unterscheidet sich je nach Wahl der Methode zur Lösung des multikriteriellen Optimierungsproblems und je nach Diskretisierung. Einzelheiten hierzu finden sich in der Dissertation1.

Methode: eval_grad_f

bool Allokation_nlp::eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f)
  vector<double>* Gradient_F=(*DMO).Gradient_Zielfunktion(x);

Der Gradient wird im Allokationsmodul auf zwei Weisen berechnet:

Methode: eval_g

bool Allokation_nlp::eval_g(Index n, const Number* x, bool new_x, Index m, Number* g)

Die Nebenbedingungen enthalten die Restriktionen der Antriebe; beispielsweise bezüglich ihrer Leistung (vgl. Dissertation).

//VSP
  for(int i=0;i<m_VSP;i++){
    for(int k=0;k<z;k++){
      (*Nebenbedingungen)[k+i*z]=(pow((*VSP_x)[k+i*z],2)+pow((*VSP_y)[k+i*z],2))/(pow(Steigung_max[i],2));    //Pitch
      (*Nebenbedingungen)[k+i*z+m_VSP*z]=((*VSP_n)[k+i*z]-Drehzahl_min_VSP[i])/(Drehzahl_max_VSP[i]-Drehzahl_min_VSP[i]);    //Drehzahl
      (*Nebenbedingungen)[k+i*z+2*m_VSP*z]=((*VSP_Moment)[k+i*z]-Moment_min_VSP[i])/(Moment_max_VSP[i]-Moment_min_VSP[i]);
      (*Nebenbedingungen)[k+i*z+3*m_VSP*z]=((*Leistung_VSP)[k+i*z]-Power_min_VSP[i])/(Power_max_VSP[i]-Power_min_VSP[i]);    //Leistung
    }
  }

//Azimuth
  for(int i=0;i<m_Azi;i++){
    for(int k=0;k<z;k++){
      (*Nebenbedingungen)[k+i*z+4*m_VSP*z]=((*Azi_n)[k+i*z]-Drehzahl_min_Azi[i])/(Drehzahl_max_Azi[i]-Drehzahl_min_Azi[i]);    //Drehzahl
      (*Nebenbedingungen)[k+i*z+4*m_VSP*z+m_Azi*z]=((*Azi_Moment)[k+i*z]-Moment_min_Azi[i])/(Moment_max_Azi[i]-Moment_min_Azi[i]);      
      (*Nebenbedingungen)[k+i*z+4*m_VSP*z+2*m_Azi*z]=((*Leistung_Azi)[k+i*z]-Power_min_Azi[i])/(Power_max_Azi[i]-Power_min_Azi[i]);   //Leistung
    }
  }

//Bugstrahler
  for(int i=0; i<m_Bug;i++){
    for (int k=0; k<z;k++){
      (*Nebenbedingungen)[k+i*z+4*m_VSP*z+3*m_Azi*z]=((*Bug_n)[k+i*z]+1)/(1+1);
      (*Nebenbedingungen)[k+i*z+4*m_VSP*z+3*m_Azi*z+m_Bug*z]=(*Leistung_Bug)[k+i*z]/Power_max_Bug[i];
    }
  }
vector<double>* Nebenbedingung_g=(*DMO).Nebenbedingung(x);

Methode: eval_jac_g

bool Allokation_nlp::eval_jac_g(Index n, const Number* x, bool new_x,
                           Index m, Index nele_jac, Index* iRow, Index *jCol,
                           Number* values)

Hier wird zunächst die Struktur der Jacobimatrix festgelegt, also die Zeile und Spalten in denen ein Eintrag ungleich Null steht.

    vector<double>* Jacobi_ROW=(*DMO).Struktur_Jacobi_ROW();  
    vector<double>* Jacobi_COL=(*DMO).Struktur_Jacobi_COL();

Aufgrund der Modellierung der Steuerung durch konstante B-Splines ergibt sich eine dünnbesetzte Jacobimatrix. Die Dreiecksstruktur der Matrix ergibt sich aus der Berechnung der Zustandsvariablen, wie beispielsweise der Drehzahl

\[n_{k+1}=n_0+h\tilde{n}\sum_{i=0}^k\omega_i \]

mit der Ableitung nach \(\omega_j\)

\[n_{k+1,\omega_j}=h\tilde{n}\chi_{j\leq k}. \]

Für die Berechnung der Jacobimatrix werden die Nebenbedingungen zunächst nach ihren Zustandsvariablen und dann nach ihren Steuervariablen abgeleitet.

Beispiel: Leistung eines Azimutpropellers

Die Nebenbedingungen, die die Leistung des Azimutpropellers zu den Zeitpunkten \(t_1\) und \(t_2\) beschränken haben die Form:

\[g(t_k)=\frac{c_{\text{Azi}}n(t_k)^3-P_{\min}}{P_{\max}-P_{\min}},\ \ k=1, 2. \]

Die Konstante \(c_{\text{Azi}}\) ist dabei die propellerspezifische Konstante für die Leistung (vgl. Dissertation1, Kapitel 3). Abgeleitet nach den Zustandsvariablen des Azimutpropellers gilt:

\[g(t_k)_{\alpha}=0,\ \ k=1, 2\\ g(t_k)_{n}=\frac{3c_{\text{Azi}}n(t_k)^2}{P_{\max}-P_{\min}},\ \ k=1, 2. \]

Weiter gilt für die Zustandsvariablen abgeleitet nach den Steuervariablen:

\[n(t_k)_{\omega_0}=1,\ \ k=1, 2\\ n(t_1)_{\omega_1}=0,\\ n(t_2)_{\omega_1}=1. \]\[\begin{equation*} \Rightarrow J_g:=\begin{pmatrix} 0&0&x&0\\ 0&0&x&x \end{pmatrix} \end{equation*} \]

Alle Nebenbedingungen sind ausführlich in der Dissertation1 in Kapitel 4 beschrieben.

Die Struktur der Jacobimatrix wird nun am Beispiel von je zwei VSP, zwei Azimutpropellern und zwei Bugstrahler dargestellt. Die roten Dreiecke entstehen, falls Interaktionen verwendet werden.

Danach werden die Einträge der Matrix berechnet.

    vector<double>* Jacobi_value=(*DMO).Jacobi_Inhalt(x);

Wie bei der Berechnung des Gradienten der Zielfunktion werden auch hier die Einträge analytisch bestimmt, falls keine Interaktionen beachtet werden; ansonsten erfolgt die Berechnung über finite Differenzen.

Methode: finalize_solution

void Allokation_nlp::finalize_solution(SolverReturn status,
                                  Index n, const Number* x, const Number* z_L, const Number* z_U,
                                  Index m, const Number* g, const Number* lambda,
                                  Number obj_value,
                                  const IpoptData* ip_data,
                                  IpoptCalculatedQuantities* ip_cq)

Für eine Folge an Allokationsauswertungen muss nach jeder Optimierung die aktuelle Steuerung zwischengespeichert werden.

(*DMO).Save_Startwerte(x);  //fuer eine zusammenhaengende Folge an Auswertungen

Im Falle der lexikographischen Optimierung ist zudem das Speichern der Lösung nach der ersten Optimierung notwendig um einen zulässigen Startwert in der zweiten Optimierung zu gewährleisten.

(*DMO).Save_Steuerung(x);    //Speichert die Steuerung x

Des Weiteren muss der Propellerschub aus der ersten Optimierung gesichert werden, um die zusätzliche Nebenbedingung in der zweiten Optimierung aufstellen zu können.

(*DMO).Save_Tau_thruster(x); //Berechnet die Kraft tau bei der Lexikographischen Optimierung

All diese Ergebnisse befinden sich im Ordner Zwischenergebnisse und werden nach der Optimierung wieder überschrieben.

Zuletzt werden noch Resultate wie Propellerschub der einzelnen Antriebe, Leistung, Moment usw. gespeichert. Diese sind hauptsächlich bei Simulationen von Interesse.

(*DMO).Save_Ergebnisse(x); 

Zusätzliche Funktionen

Die für die Zielfunktion und die Nebenbedingungen notwendigen Propellerschübe und Momente bzw. Leistung,` werden mit den Funktionen

Thrust_VSP(vector<double> pitch_x, vector<double> pitch_y, vector<double> n_VSP)
Moment_VSP(vector<double> pitch_x, vector<double> pitch_y, vector<double> n_VSP)
Power_VSP(vector<double> pitch_x, vector<double> pitch_y, vector<double> n_VSP)
Thrust_Azimuth(vector<double> alpha, vector<double> n_Azi)
Moment_Azi(vector<double> n_Azi)
Power_Azi(vector<double> n_Azi)
Thrust_Bug(vector<double> n_Bug)
Power_Bug(vector<double> n_Bug)

berechnet. Dabei unterscheiden sich die Funktionen im Falle einer Interaktion von den Funktionen in denen keine Interaktion beachtet wird. Die genaue Modellierung der Funktionen findet sich in der Dissertation1 in Kapitel 4.

Variante des Allokationsmoduls

Das Ziel in der Allokation ist die Verteilung einer geforderten Kraft auf die Antriebe, sodass die geforderte Kraft möglichst gut und schnell und zudem leistungseffizient erzeugt wird. Dabei kann leistungseffizient bedeuten, dass die Gesamtleistung der Antriebe minimal werden soll, oder dass jeder einzelne Propeller möglichst wenig Leistung benötigt. Um beiden Interpretationen gerecht zu werden, gibt es neben dem Ordner Allokationsmodul auch Allokationsmodul_Variante. Führt man die Datei Allokation_Ipopt im Ordner Allokationsmodul aus, so wird die Leistung jedes einzelnen Propellers minimiert, im Ordner “Allokationsmodul_Variante” dagegen kommt es zu einer Mimimierung der Gesamtleistung.

Footnotes