WICHTIG: Der Betrieb von goMatlab.de wird privat finanziert fortgesetzt. - Mehr Infos...

Mein MATLAB Forum - goMatlab.de

Mein MATLAB Forum

  
Gast > Registrieren       Autologin?   

Partner:




Forum
      Option
[Erweitert]
  • Diese Seite per Mail weiterempfehlen
     


Gehe zu:  
Neues Thema eröffnen Neue Antwort erstellen

Ballistisches Problem
 

schalom
Forum-Newbie

Forum-Newbie

Beiträge: 8
Anmeldedatum: 10.01.12
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 14.01.2012, 15:39     Titel: Ballistisches Problem
  Antworten mit Zitat      
Hallo Leute,

ich muss mit Matlab die Flugbahn eines Geschosses unter Berücksichtigung des Luftwiderstands darstellen.

Die Differentialgleichung sieht so aus:

mx''=-10^-3*sqrt(x'^2+y'^2)*x'
my''=-10^-3*sqrt(x'^2+y'^2)*y'-10

ich hab folgenden Ansatz mit Matlab versucht:

Code:
%x1=x, x2=x', x1'=x2, x2'=10^-3*sqrt(x'^2+y'^2)*x'
%x3=y, x4=y', x3'=x4, x4'=10^-3*sqrt(x'^2+y'^2)*y'-mg



function dydx = aufgabe1(x,y)
dydx=[x(4);-10*10^-3*(sqrt(x(1)^2+x(3)^2))*x(3)-10,x(2);-10*10^-3*(sqrt((x(1))^2+(x(3))^2))*x(1)];

 


In einer 2 m-file habe ich folgendes reingeschrieben:


Code:
 
x0=0;    % Startpunkt
y0=1;    % Startwert
xend=10;  % Endpunkt
[x,y] = ode45('aufgabe1', [x0, xend], y0);
% Spaltenvektor x enthaelt die Stuetzstellen,
% Spaltenvektor y die dort vorliegenden Werte der Integralkurve.
lenx=length(x);
plot(x,y,x,zeros(1,lenx),'o')
%
figure
ode45(@aufgabe1, [x0, xend], y0)



Leider funktioniert das nicht. Folgende Fehlermeldung erscheint:

Code:
??? Attempted to access x(4); index out of bounds because numel(x)=1.

Error in ==> aufgabe1 at 11
dydx=[x(4);-10*10^-3*(sqrt(x(1)^2+x(3)^2))*x(3)-10,x(2);-10*10^-3*(sqrt((x(1))^2+(x(3))^2))*x(1)];

Error in ==> odearguments at 109
f0 = feval(ode,t0,y0,args{:});   % ODE15I sets args{1} to yp0.

Error in ==> ode45 at 173
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...

Error in ==> aufgabe at 11
[x,y] = ode45('aufgabe1', [x0, xend], y0);
 



Kann mir jemand sagen, wo der Fehler liegt?


Vielen Dank
Private Nachricht senden Benutzer-Profile anzeigen


Harald
Forum-Meister

Forum-Meister

Beiträge: 24.492
Anmeldedatum: 26.03.09
Wohnort: Nähe München
Version: ab 2017b
     Beitrag Verfasst am: 14.01.2012, 15:59     Titel:
  Antworten mit Zitat      
Hallo,

deine DGL hat 4 Dimensionen, du musst also als Startvektor (y0) auch einen Vektor mit 4 Elementen vorgeben.

Bei deiner Funktion solltest du auf , und ; aufpassen, und auf die Reihenfolge der Argumente achten, die ode45 vorgibt. Du müsstest immer y(1) bis y(4) statt x(1) bis x(4) verwenden. Mir scheint zudem, dass du 2 und 4 statt 1 und 3 verwenden musst, um deine DGL korrekt umzusetzen. Jedenfalls stimmt dein Ansatz nicht mit deiner Implementierung überein.

Grüße,
Harald
Private Nachricht senden Benutzer-Profile anzeigen
 
schalom
Themenstarter

 Forum-Newbie

Forum-Newbie

Beiträge: 8
Anmeldedatum: 10.01.12
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 14.01.2012, 16:25     Titel:
  Antworten mit Zitat      
Dank für deine Hilfe, aber es funktioniert immer noch nicht richtig Sad

Die Funktion habe ich jetzt geändert:

Code:

function dydx =

function dydx = aufgabe1(x,y)
dydx=[y(4),-10*10^-3*(sqrt(y(2)^2+y(4)^2))*y(4)-10;y(2),-10*10^-3*(sqrt((y(2))^2+(y(4))^2))*y(2)];
 


Die m-File sieht nun so aus:

Code:
 
 
x0=0;    % Startpunkt
y0=[0,0,0,0] % Startwert
yend=[10,10,10,10];  
  % Endpunkt
[x,y] = ode45('aufgabe1', [y0, yend], [x0]);

leny=length(y);
plot(x,y,x,zeros(1,lenx),'o')
%
figure
ode45(@aufgabe1, [y0, yend], [x0])

 



Als Fehlermeldung erscheint nun folgendes:
Code:
??? Error using ==> odearguments at 105
The entries in tspan must strictly increase or decrease.

Error in ==> ode45 at 173
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...

Error in ==> aufgabe at 6
[x,y] = ode45('aufgabe1', [y0, yend], [x0]);
 



Was kann ich verbessern?
Private Nachricht senden Benutzer-Profile anzeigen
 
Harald
Forum-Meister

Forum-Meister

Beiträge: 24.492
Anmeldedatum: 26.03.09
Wohnort: Nähe München
Version: ab 2017b
     Beitrag Verfasst am: 14.01.2012, 18:50     Titel:
  Antworten mit Zitat      
Hallo,

du kannst aufhören, x und y zu verwechseln Wink

Als Ende musst du das Ende der Simulationszeit angeben, keinen Endzustand. Wenn du einen Endzustand angeben möchtest, musst du ein anderes Verfahren anwenden, z.B.
Code:


Grüße,
Harald
Private Nachricht senden Benutzer-Profile anzeigen
 
schalom
Themenstarter

 Forum-Newbie

Forum-Newbie

Beiträge: 8
Anmeldedatum: 10.01.12
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 15.01.2012, 17:04     Titel:
  Antworten mit Zitat      
Also ich bin nach dem Beispiel aus der Matlab Help vorgegangen und hab nun folgendes:

Code:
function dy = rigid(x,y)
dy = zeros(2,1);    % a column vector
dy(2) = -10*10^-3*(sqrt((y(2))^2+(y(4))^2))*y(2);
dy(4) = -10*10^-3*(sqrt(y(2)^2+y(4)^2))*y(4)-10;
 



In einer zweiten m-file hab ich folgendes reingeschrieben:


Code:
options = odeset('RelTol',1e-4,'AbsTol',[1e-3 1e-3 1e-3 1e-4]);
[x,Y] = ode45(@rigid,[0 10],[0 0 0 0],options);

lenx=length(Y);
plot(x,Y,x,zeros(1,lenx),'o')

figure
ode45(@rigid,[0 10],[10 0 0 0],options)



Fehlermeldungen bekomme ich keine mehr. Auch erscheint jetzt ein Plot. Leider sieht dieser Plot nicht wie erwartet aus. Die geplottete Kurve müsste doch mehr oder weniger (wegen dem Luftwiderstand) parabelförmig aussehen. Tut sie aber bei mir nicht.

Kann mir einer sagen wo der Fehler liegt?

Danke Very Happy
Private Nachricht senden Benutzer-Profile anzeigen
 
Jan S
Moderator

Moderator

Beiträge: 11.057
Anmeldedatum: 08.07.10
Wohnort: Heidelberg
Version: 2009a, 2016b
     Beitrag Verfasst am: 15.01.2012, 17:22     Titel:
  Antworten mit Zitat      
Hallo schalom,
In Deiner rigid Funktion möchtest Du doch einen Vektor nmit 4 Elementen erzeugen, oder? Dann ist "zeros(2,1)" keine gute Wahl. Sollen die Elemente dy(1) und dy(3) wirklich 0 sein?!

Bemerkung: -10*10^-3 benötigt ein paar teure Operationen, -10e-3 dagegen keine.

Gruß, Jan

Zuletzt bearbeitet von Jan S am 15.01.2012, 23:04, insgesamt einmal bearbeitet
Private Nachricht senden Benutzer-Profile anzeigen
 
schalom
Themenstarter

 Forum-Newbie

Forum-Newbie

Beiträge: 8
Anmeldedatum: 10.01.12
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 15.01.2012, 17:28     Titel:
  Antworten mit Zitat      
Hallo,

also ich hab keine konkreten Werte vorgegeben.

dy1 wäre ja die Geschwindigkeit in x-Richtung und die wäre ja in dem Moment in dem ich das Objekt werfe null. dy3 wäre die Anfangsgeschwindigkeit in y-Richtung und die wäre ebenso null oder?

Aber auch wenn ich die Anfangswerte ändere, ändert sich am Verlauf der Funktion sehr wenig.
Private Nachricht senden Benutzer-Profile anzeigen
 
Harald
Forum-Meister

Forum-Meister

Beiträge: 24.492
Anmeldedatum: 26.03.09
Wohnort: Nähe München
Version: ab 2017b
     Beitrag Verfasst am: 15.01.2012, 17:30     Titel:
  Antworten mit Zitat      
Hallo,

du rechnest in beiden Fällen mit Anfangsgeschwindigkeit 0, also quasi einem senkrechten Fall.

Zudem hast du in deiner letzten Variante plötzlich 2 DGLen weniger.

So siehts gut aus:
Code:
[x,Y] = ode45(@aufgabe1,[0 10],[0 10 0 0],options);

plot(Y(:,1), Y(:,3))


mit
Code:
function dydx = aufgabe1(x,y)
dydx=[y(2);
    -10*10^-3*(sqrt(y(2)^2+y(4)^2))*y(2);
    y(4);
    -10*10^-3*(sqrt(y(2)^2+y(4)^2))*y(4)-10];


Grüße,
Harald
Private Nachricht senden Benutzer-Profile anzeigen
 
schalom
Themenstarter

 Forum-Newbie

Forum-Newbie

Beiträge: 8
Anmeldedatum: 10.01.12
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 15.01.2012, 17:43     Titel:
  Antworten mit Zitat      
Vielen Dank Very Happy
Private Nachricht senden Benutzer-Profile anzeigen
 
Neues Thema eröffnen Neue Antwort erstellen



Einstellungen und Berechtigungen
Beiträge der letzten Zeit anzeigen:

 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
.





 Impressum  | Nutzungsbedingungen  | Datenschutz | FAQ | goMatlab RSS Button RSS

Hosted by:


Copyright © 2007 - 2024 goMatlab.de | Dies ist keine offizielle Website der Firma The Mathworks

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.