Verfasst am: 18.02.2016, 12:29
Titel: Mit ode45 DGL 1.Ordnung in Matrixform lösen
Hallo nochmal!
Ich probiere es noch einmal in einem neuen Thread, da ich kurz vorm Verzweifeln bin
Ich möchte gerne eine DGL 1. Ordnung, die Matritzen unterschiedlicher Dimensionen enthält, mithilfe des Befehls ode45 lösen. Die Differentialgleichung lautet wie folgt:
dy/dt=(inv(L)*A*E(t))' + inv(L)*D*(y(t))'
Informationen zu den Matritzen:
1. dy/dt ist ein 4x1 Vektor, wobei für jeden dieser 4 Einträge ydot1, ydot2, ydot3 und ydot4 die DGL gelöst werden soll. Als Ergebnis soll dann eine nt x4 Matrix rauskommen, sodass die erste Zeile beispielsweise die Ergebnisse aller ydot(1-4) zum Zeitounkt=0 enthält.
2. L ist eine 4x4 Matrix und konstant
3. A ist eine 4x6 Matrix und konstant
4. D ist eine 4x4 Matrix und konstant
5. E(t) ist eine 6x nt Matrix, wobei nt für die Zeitschritte steht
Ich habe nun mehrere Fragen:
1. Die Matrix E ist in sofern zeitabhängig, dass jede Spalte dieses Vektors einem Zeitpunkt zugeordnet ist. Deswegen hat die Matrix auch nt Spalten. Wie baue ich diese Information in die ode45 ein? Es müsste also für jeden Zeitschritt die entsprechende Spalte verrechnet werden.
2. Ich habe zudem gelesen, dass für Y ein Spaltenvektor rauskommen muss, wobei jede Zeile einem Zeitschritt zugeordnet ist. Wenn ich die entsprechenden Matritzen transponiere, wie oben, dann bekomme ich für die Ergebnis-Matrix die Form nt x4 (ohne Transponieren, wäre die Matrix um 90° gedreht und dann würde nicht jede Zeile einem Zeitschritt entsprechen, wie es gefordert ist). Wie kann ich diese Ergebnis-Matrix vorher so als Array definieren, dass einzelne Vektoren abgespeichert werden, aber eben vier Stück für ydot(1-4)?
3: Überhaupt habe ich das Problem mit den Dimensionen. Welche Dimension wird für y(t) angenommen und kann ich das steuern? Sie müsste in irgendeiner Form die Dimensionen nt x4 bzw. nach dem Transponieren 4x nt haben.
4. Welche Dimension hat dann die Anfangsbedingung y0? wäre dann y0=[0,0,0,0] für die 4 zu berechnenden Variablen?
Hier nocheinmal die Codes. Zunächst das Hauptskript:
dt = 1e-5; %time stepsize
tstart = 0; %starting time
tend = 1e-1; %total time
tr = 1e-3; %loading time (ramp)
te =[tstart:dt:tend]';
nt = length(te); %number time steps
was soll denn dieses E sein? ist das eine zeitabhängigkeit? wie man die benutzt steht doch in der doc unter example 3 ich würd mich daran halten. oder soll das was anderes sein? E sollte also 6x1 sein und y0 4x1. irgendwie finde ich auch komisch das 5 einträge von E dann 0 sind. da kann man bestimmt irgendwo auch was einspaaren dann.
_________________
Wenn man E plottet zeigt es eine zunächst linear ansteigende Rampe. Ab einem bestimmten Zeitpunkt ist E dann konstant. Dass die anderen Einträge gleich Null sind ist nicht immer der Fall, deswegen kann E auch in keinem Fall vereinfacht werden.
Und wie bereits gesagt: Die 6 Einträge pro Spalte zeigen einen Zustand zu einem bestimmten Zeitschritt nt.
Und ich habe mir schon unzähliges dazu durchgelsen, aber ich verstehe es nicht. Deswegen stelle ich die Frage hier ja.
dt = 1e-5; %time stepsize
tstart = 0; %starting time
tend = 1e-1; %total time
tr = 1e-3; %loading time (ramp)
te =[tstart:dt:tend];
nt = length(te); %number time steps
Ok, es ist schonmal super, dass es durchläuft, also vielen herzlichen Dank dafür!
Eine Frage habe ich aber noch: ist es richtig E zu interpolieren bzw. gibt das nicht ein falsches Ergebnis? E liefert nur bis zum Teitpunkt tr linear ansteigende Werte, ab dann beträgt der Wert immer nur 0.01. Im Prinzip sind das ja dann zwei Funktionen, wenn man so will. Kann E dann als eine einzige Funktion interpoliert werden?
Kann ich das irgendwie mit einer if-Schleife lösen?
Beispielsweise so in der Art:
Code:
function dydt = myode3(tr,t,y,E,A,D,L,te) if t<tr
E=interp1(te,E',t)';
else
E(1,:)=0.01;
end
dydt = (inv(L)*A*E)+inv(L)*D*(y); % Evaluate ODE at time t end
: ist es richtig E zu interpolieren bzw. gibt das nicht ein falsches Ergebnis? E liefert nur bis zum Teitpunkt tr linear ansteigende
solange der anstieg nicht unendlich ist und ist es doch egal für die interpolation.
du kannt auch eine funktion da hin schreiben für E
_________________
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.