Ich möchte nun mein System erweitern und die Masse zwischen 2 Stabfedern einspannen, die jeweils nur Druckkräfte auf die Masse ausüben. Die zusätzlichen Kräfte durch die Stabfedern sind Fl(x) und Fr(x)
Zur Berechnung dieser Kräfte habe ich ein Werkstoffmodell zur Verfügung, welches mir in Abhängigkeit von der Dehnung eps der Stabfedern eine Spannung sigma ausgibt.
Der Spannung-Dehnungs-Verlauf dieser Stabfedern beschreibt dabei eine Hysterese.
Über den Zusammenhang Kraft = Spannung*Fläche möchte ich dann Fr und Fl berechnen.
Die Spannung der jeweiligen Stabfeder wird durch ein zeitinkrementelles Prädiktor-Korrektor-Verfahren berechnet. Bei diesem Verfahren wird die Spannung zum Zeitpunkt tn aus einem Prädiktor und einem Korrektor-Schritt berechnet.
Die Herausforderung liegt nun in der Umsetzung dieses Vorgehens. Wie verknüpfe ich das Lösungsverfahren der DGL (ode45) im Programm "solver" mit dem Prädiktor-Korrektor-Verfahren im Programm "Werkstoffmodell"?
die Klasse Werkstoffmodell ist sehr umfangreich. Ich habe keine Ahnung, wie sie genutzt werden soll. Daher kann ich erst recht nicht sagen, wie sie mit ode45 verknüpft werden kann.
Solltest du den Autor der Klasse kennen, spreche doch mal ihn darauf an.
die Vorgehensweise des Programms soll wie folgt aussehen:
Zu Beginn wird die Funktion "Werkstoffmodell.fnSetPrestress" ausgeführt. Diese berechnet die Vorspannung und -dehnung in den Stabfedern zum Zeitpunkt t=0. Im Anschluss daran wird die Verlagerung x der Masse zum Zeitpunkt t=0 berechnet. In meinem Fall beträgt diese 0, da noch keine äußere Kräften wirken und die Startbedingung x(0)=0 ist.
Im nächsten Zeitschritt (t1) wirkt dann die Anregungskraft F(t1) auf das System, wodurch es ausgelenkt wird. Die Auslenkung der Masse führt dazu, dass die rechte Stabfeder weiter gestaucht und die linke entlastet wird. Dies führt dazu, dass die Kraft Fr(t1)>Fl(t1) ist.
Zur Berechnung von Fr(t1) und Fl(t1) wird zunächst die Dehnung eps_r(t1) und eps_l(t1) berechnet.
Die Funktion "Werkstoffmodell.fnCurrentState(eps_r_n,eps_l_n,t_n)" berechnet dann ein delta_eps und delta_t aus der Differenz der jeweiligen Werte zum aktuellen Zeitpunkt t1 und dem vorherigen Zeitpunkt t0. Diese Werte werden dann in die Funktionen "Werkstoffmodell.fnRunPredictorStep; Werkstoffmodell.fnEvaluateFlowCondition; Werkstoffmodell.fnUpdateTemperature;" eingesetzt, um die Spannung der Stabfedern sig_r und sig_l zum Zeitpunkt t1 zu berechnen.
Die Stabfederkräfte Fr und Fl können nun über F=sig*A berechnet werden und das Runge-Kutta Verfahren kann bis zum Zeitpunkt t1 die DGL lösen und die Verlagerung der Masse x(t1) ausgeben.
Im letzten Schritt wird die Funktion "Werkstoffmodell.UpdateStateVariables" aufgerufen, um die Zustandsgrößen abzuspeichern.
Nachdem die zuvor beschriebenen Funktionen ausgeführt wurden, wird der ganze Ablauf für den nächsten Zeitschritt (t2) durchgeführt, solange bis t_end erreicht wird.
Ich hoffe das Vorgehen ist so wie beschrieben nachvollziehbar. Über eine Hilfestellung zur Implementierung dieses Vorgehens würde ich mich sehr freuen.
Hallo,
Ich möchte das Thema noch einmal pushen, da ich bisher keine Lösung zu meinem Problem gefunden habe.
Ist das von mir beschriebene Vorgehen so überhaupt mit dem ode-solver (Runge-Kutta) möglich?
An der Umsetzung scheitert es aktuell daran, dass die Klasse "Werkstoffmodell" zu jedem Zeitpunkt aufgerufen werden muss, wohingegen die Funktion Ode45 einmalig für den gesamten Zeitbereich aufgerufen wird. Liege ich damit richtig?
Um dieses Problem zu entgehen, könnte ich alternativ zum ode-solver die Newmark-beta-Methode implementieren. Diese Methode wird wie die Klasse "Werkstoffmodell" zu jedem Zeitpunkt aufgerufen. Macht diese Methode in meinem Fall mehr Sinn?
Über eine Rückmeldung zu diesem Thema würde ich mich sehr freuen.
An der Umsetzung scheitert es aktuell daran, dass die Klasse "Werkstoffmodell" zu jedem Zeitpunkt aufgerufen werden muss, wohingegen die Funktion Ode45 einmalig für den gesamten Zeitbereich aufgerufen wird. Liege ich damit richtig?
Du kannst die Klasse ja in der "rechten Seite" der DGL aufrufen, und dann wird sie auch in jedem Zeitschritt aufgerufen. Problematisch kann sein, dass bei ode45 auch mal ein Schritt abgelehnt wird, wenn Fehlertoleranzen nicht eingehalten werden.
Erster Schritt dürfte sein, die "rechte Seite" mal komplett in eine Funktion auszulagern statt ein Anonymous Function Handle zu verwenden.
1.) Wie du bereits erwähnt hast rechnet ode45 je nach Fehlertoleranz mit unterschiedlichen Zeitschritten dt. Kann ich diese Zeitschritte auch als konstant vorgeben?
Dieser Ansatz funktioniert nicht so ganz, da dt von ode45 verändert wird:
Code:
t=0:obj.dt:obj.t_end; % time scale [obj.t1,obj.y1] = ode45(@(t,x)[x(2); (-x(2)*obj.c-obj.k*x(1)+Fresfun(x,t,obj.Werkstoffmodell,obj)+obj.F(t))/obj.m] ,t [00]); % Neue DGL
2.) Wie kann ich die von ode45 verwendeten Zeitschritte dt abspeichern bzw. auslesen? Für meine Klasse "Werkstoffmodell" benötige ich die genaue Schrittweite.
mein Fehler, es gibt nur ein MaxStep.
Das wäre aber so oder so nur eine Notlösung gewesen: es ist ja nicht sinnvoll, einen Löser mit variabler Schrittweite dazu zu zwingen, eine konstante Schrittweite beizubehalten.
Während des Testens dieser Verfahren musste ich allerdings feststellen, dass ich mit diesen solvern nur sinnvolle Ergebnisse bekomme, wenn die Masse mind. 100 kg beträgt. Für meinen Anwendungsfall möchte ich allerdings eine Masse von 0.004 kg vorgeben. In diesem Fall divergiert allerdings die Simulation für alle fixed step solver und die Ergebnisse sind nicht mehr zu gebrauchen. Auch eine Reduzierung des Zeitschritts von 1/1000 s auf 1/100000 s konnte dieses Problem nicht beheben.
Daraufhin habe ich wieder mit variable-step solver herumprobiert und konnte dabei gute Ergebnisse bei der Verwendung der stiff-ode-solver (ode15s,ode23s,ode23tb) generieren. Die nonstiff-ode-solver (ode45,ode23,ode113) konnte ich ebenfalls nicht verwenden, da die Simulationsdauer sehr hoch war.
Da es sich in diesem Fall offensichtlich um ein numerisches Problem handelt, bin ich um jeder weitere Hilfe sehr dankbar.
Mir fehlt es momentan an der systematischen Herangehensweise an dieses Problem. Handelt es sich beispielsweise bei meinem System um eine steife DGL? Im Idealfall würde ich gerne wieder auf ein fixed-step-solver umschwenken. Welche Möglichkeiten (z.B. andere numerische Verfahren) bieten sich dafür an?
Grüße
Holger Gerach
Gast
Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
Verfasst am: 17.05.2017, 23:36
Titel:
Hallo,
Ich hätte einen Extended Kalman Filter verwendet, den man leicht in einer m-function implementieren kann. Bei der gerigen Systemordnung dürfte das gut gehen. Alternative wäre ein Unscented Kalman Filter, falls die Schätzung nicht erwartungstreu genug ist.
Nun weiß ich allerdings nicht, wie ich diesen Filter in mein Programm implementieren kann. Habe ich es richtig verstanden, dass ich den Kalman-Filter verwende, nachdem ich meine DGL mit einem ode-solver gelöst habe?
Es wäre sehr hilfreich, wenn du mir eine Hilfestellung geben könntest.
hat noch jemand einen guten Hinweis zu meiner Problemstellung?
In der Literatur wird das von mir beschriebene System (nonlinear dynamics of a sdof oscillator with hysteresis) häufig durch einen "Newmark implicit algorithm" in Verbindung mit der "Newton-Raphson-Iteration" gelöst. Den Algorithmus habe ich mal in den Anhang gepackt.
Liefert dieses Verfahren überhaupt Vorteile gegenüber den in Matlab bereits implementiertetn ode-solvern (Runge-Kutta)?
Du kannst Beiträge in dieses Forum schreiben. Du kannst auf Beiträge in diesem Forum antworten. Du kannst deine Beiträge in diesem Forum nicht bearbeiten. Du kannst deine Beiträge in diesem Forum nicht löschen. Du kannst an Umfragen in diesem Forum nicht mitmachen. Du kannst Dateien in diesem Forum posten Du kannst Dateien in diesem Forum herunterladen
MATLAB, Simulink, Stateflow, Handle Graphics, Real-Time Workshop, SimBiology, SimHydraulics, SimEvents, and xPC TargetBox are registered trademarks and The MathWorks, the L-shaped membrane logo, and Embedded MATLAB are trademarks of The MathWorks, Inc.