Subsysteme: Ereignissteuerung

Beispiel 2:

Beispiel 3:

   y1'(t) = -y1(t)*(alpha*y2 + beta) + gamma
   y2'(t) = y2*(rho*y1 - sigma) + tau*(1+y1)

   y1(0) = -1,  y2(0) = 0

   alpha=1.5e-18,   beta=2.5e-6,   gamma=2.1e-6, 
   rho=0.6,         sigma=0.18,    tau=0.016

Ein Simulink-Modell zu diesem Anfangswertproblem ist: 

• Relative tolerance = 1.0e-6
• Absolute tolerance = 1.0e-8

• die Schrittweiten,
• die Anzahl der Schritte,
• den Quotienten der Eigenwerte der Jacobi-Matrix,
• die Durchlaufzeiten

Zusätzlich umfasst das Modell ein getriggertes Subsystem, das anzeigt:

• die Durchlaufzeiten (Echtzeit)
• den Quotienten von betragsgrößtem und betragskleinstem Eigenwert der Jacobi-Matrix. Diesr Quotient wird nur berechnet, solange die Realteile beider Eigenwerte negativ sind. Der Quotient wird 'Steifigkeit' oder auch 'Steifigkeitsquotient' des Systems genannt.

• Amplitude=1
• Period=2
• Pulse width=1
• Phase delay=0
• Sample time=0.5e4

rubytime liefert die Durchlaufszeit in Sekunden. Matlab kennt die beiden Kommandos tic und toc. Ersteres initialisiert eine interne Zeitmessung, letzteres liefert die Zeitdifferenz (Echtzeit) zwischen Zeitpunkt des toc-Kommandos und dem vorangegangenen tic-Kommando. Um diese Zeitmessung zu ermöglichen, wird:

• das Subsystem mit dem Initialisierungskommando 'tic;' gestartet
• umfasst die M-Function rubytime das Kommando 'toc;'

     function t = rubytime(u)
     t = toc;

eigjac liefert den Steifigkeitsquotienten des Systems und wird nur berechnet, wenn die Realteile beider Eigenwerte kleiner als 0 sind. Die M-Function lautet:

   function q = eigjac(u)
   % Jacobi-Matrix
   J = [ -(1.5e-18*u(2)+2.5e-6), -1.5e-18*u(1); 
          0.6*u(2)+0.016,       0.6*u(1)-0.18 ];
   % Eigenwerte(J)
   w = eig(J);
   % wenn alle Realteile < 0
   if ( real(w) < 0 )
      % Verhaeltnis der Eigenwerte
      q = real(w(1))/real(w(2));
      if ( q < 1 )
         q = 1/q;
      end
   else
      % in diesem Falle ist der Quotient bedeutungslos
      q = 0;
   end

Ergebnisse des Probelaufes:
Wir wählen als Simulationsparameter für den Probelauf:

• Relative tolerance = 1.0e-6
• Absolute tolerance = 1.0e-8
• Maximale Schrittweite = 0.5e4 (vorgegeben durch die sample time der diskreten Blöcke)
• Simulationsintervall von 0 bis 4.0e4

• die Werte des Steifigkeitsquotienten, die etwa bei 3e5 liegen,

  

• die Durchlauszeiten bei unterschiedlichen Integrationsverfahren

  Mit dem Adams-Integrator ode113 dauert die Simulation auf einem 200 MHz Pentium-Rechner etwa 10 Sekunden:

  

  Die Durchlaufszeit beim Dormand-Prince-Verfahren liegt bei etwa 8 Sekunden.
  Deutlich darunter liegen die Durchlaufszeiten bei steifen Integratoren. Das ODE15s-Verfahren benötigt etwa 0.6 Sekunden, die zu Beginn der Simulation anfallen:

  

  Das modifizierte Rosenbrock-Verfahren benötigt für diese Simulation etwa 3 Sekunde, die wie beim ode15s-Verfahren, im wesentlichen zu Beginn der Simulation anfallen.

• die Anzahl der Integrationsschritte bei unterschiedlichen Integrationsverfahren

  Der Adams-Integrator ode113 benötigt etwa 18000 Schritte, wobei diese über das ganze Simulationsintervall etwa gleichverteilt anfallen.

  

  Das Dormand-Prince Verfahren benötigt etwa 8000 Schritte, die über das gesamte Intervall etwa gleichverteilt anfallen.
  Der steife Integrator ode15s benötigt demgegenüber nur knapp 50 Iterationsschritte. Auch bei den Iterationsschritten fällt ein großer Teil zu Beginn der Simulation an, auch wenn das Verhältnis hier nicht so deutlich ist wie bei der Simulationsdauer.

  

  Das modifizierte Rosenbrock-Verfahren benötigt für das Problem etwa 300 Schritte, die etwa gleichverteilt über dem gesamten Simulationsintervall anfallen. Es ähnelt hier eher den nicht-steifen Integratoren als ode15s.
• die Schrittweiten der unterschiedlichen Integratoren

  Der Adams-Integrator ode113 benutzt im gesamten Intervall eine Schrittweite von etwa 2, wobei die Schrittweite mit der Dauer der Simulation leicht wächst.

  

  Das Dormand-Prince-Verfahren benutzt eine mittlere Schrittweite von 5, die mit der Simulationsdauer deutlicher wächst als beim Adams-Verfahren. Die Einbrüche bei 1, 2.5 und 3.5 sind dadurch verursacht, dass diese Ereigniszeitpunkte durch die diskreten Blöcke vorgegeben werden und schlecht in das Schrittmuster des Dormand-Prince-Verfahrens passen.

  

  Der steife Integrator ode15s benutzt eine mittlere Schrittweite von etwa 1800:

  

  Der Einbruch bei 0.5 und 3.5 ist wieder auf die Vorgabe von Ereigniszeitpunkten durch diskrete Blöcke zurückzuführen, die ganz schlecht in das Schrittweitenmuster des ode15s-Verfahrens passen. Man sieht, dass das ode15s-Verfahren geraume Zeit braucht um sich von solchen Störungen zu erholen.

  Das modifizierte Rosenbrock-Verfahren benutzt etwa Schritte der Länge 150, die im Verlaufe der Simulation erkennbar kürzer werden. Die vorgegebenen Ereigniszeitpunkte stören das Schrittweitenmuster dieses Verfahrens deutlich.

  

Wir wollen die Laufzeiten und den Steifigkeitsquotienten auch bei dieser langen Simulation beobachten, um den Aufwand dafür zu begrenzen, wollen wir diese Beobachtungen aber nur alle 0.5e6 Zeiteinheiten vornehmen. Die maximale Schrittweite passen wir dementsprechend an.

Die Integration wird in etwa 11 Sekunden durchgeführt, wobei der größte Aufwand im Intervall von 0.5e6 bis 1.0e6 getrieben werden muss.

Wie man an diesem Beispiel sieht, gilt:

• die Steifigkeit eines Problems kann sich im Verlaufe der (simulierten) Zeit ändern, dementsprechend können in unterschiedlichen Bereichen unterschiedliche Integratoren zur Behandlung des Problems angemessen sein, was hier allerdings nicht ausprobiert wurde,
• die Lösungskurven steifer Probleme sind meist nicht sonderlich auffällig. Dort wo y(2) zu schwingen anfängt, ist das Problem nicht mehr steif.