ich habe folgendes Problem: Ich möchte einen Fallschirmsprung in Matlab simulieren. Die entsprechenden Pararmeter gebe ich vorher ein wie etwa Höhe, Durchmesser Fallschirm, CW Wert etc.
Mit folgender Programmierung soll dies geschehen:
function [ dydt ] = fallschirmdgl( t,y )
global g rho cw r m;
dydt=zeros(2,1);
dydt=[y(2); g-(rho*cw*pi*r^2)/(2*m).*y(2).^2];
clear all; clc;
global g rho cw m r;
g = 9.81;
rho = 1.2;
cw = 1.4;
m = 100;
r = 5;
y0=[100;50];
tspan=0:0.01:1;
[t,y]=ode45(@fallschirmdgl,tspan,y0);
v=y(:,2);
s=y(:,1)-100;
subplot(121);
plot(t,s);
xlabel('t');
ylabel('s');
subplot(122);
plot(t,v);
xlabel('t');
ylabel('v');
Im Prinzip besteht der Fallschirmsprung aus 2 gleichen Differentialgleichungen: Die für den freien Fall von m bis Höhe X und die für den gebremsten Fall am Fallschirm von x-0m
(Typ wie oben beschrieben). Der einzigste Unterschied ist ja CW Wert.
Nun meine Frage: Wie verknüpfe ich die eine Dgl mit der anderen?
du hast ja schon Vektoren für die Funktion und die Ableitung. Du musst lediglich einen großen Vektor mit allen Funktionen und Ableitungen aufstellen, dann etwa:
Y' = f(t, Y) mit
Y(1) = x
Y(2) = x'
Y(3) = y
Y(4) = y'
ist mir noch nicht ganz klar.... aus der Hilfe für Ode45 werd ich in Bezug auf mein Beispiel noch nicht schlau daraus! Könntest du es irgendwie exakter aufzeigen? In einem lauffähigen Programm etwa? Wäre dir super Dankbar!
function[ dYdt ] = fallschirmdgl( t,Y ) global g rho cw r m;
dydt=[Y(2); g-(rho*cw*pi*r^2)/(2*m).*Y(2).^2];
dxdt = [Y(4); Formel mit Y(4)];
dYdt = [dydt; dxdt];
danke schon mal für diesen ansatz, ich hab versucht den gesamten ablauf in matlab mit deiner modifikation zu programmieren...er gibt mir allerdings einige fehler aus
ich wäre dir sehr dankbar, wenn du mir evtl einen komplettierten file zum vollständigen fallschirmsprung schicken würdest, weil ich mich nicht in der lage sehe dieses Problem zu lösen
was hältst du von dem umgekehrten Vorschlag? Du schickst, was du bist jetzt hast, und dann kann man weitersehen.
Ich werde in den nächsten Tagen wenig Zeit haben. Um das zusammenzuschreiben, würde ich aber zumindest die Anfangsbed. und die DGL brauchen, die du für die x-Richtung modellieren willst.
Ansonsten: wenn man's selber macht, lernt man mehr
%%%%%%%%%%%%%%%%%%%%%%
% Data of the jump
%%%%%%%%%%%%%%%%%%%%%%
H0 = 3000; % height of plane meter
H1 = 1500; % opening height parachut meter
Am = 0.5; % area man meter^2
Ap = 30; % area parachut meter^2
A = Am+Ap; % whole area
m = 85; % weight parachutist in kg
cw = 1.3; % cw-value
rho = 1.2; % air density kg/m^3
g = 9.81; % acceleration gravity
a1 = (cw*Am*rho)/(2*m); % air resitance before p. opens
a2 = (cw*A*rho)/(2*m); % air resitance aftzer p. opens
Es soll nun folgendes geschehen: Fallschirmspringer springt ab, es gilt Wert a1! Sobald der Fallschirmspringer den Schirm öffnet, soll die gleiche DGL abgearbetet werden, nun mit wert a2! Am Ende soll man in einem Plot die gesamte Sprungszene sehen! Danke schon mal!
leider ist mir nicht klar, wie die verschiedenen Code-Schnipsel zusammenzusetzen sind. Vielleicht mal die Fehlermeldung posten oder selbst mit dem Debugger suchen?
Du musst ja zuerst die DGL für den freien Fall lösen bis zum Öffnen bei "t_open", den du aber nicht kennst (du gibst ja die Höhe beim Öffnen vor). Die so gewonnen Werte für "x" und "v" bei "t_open" dienen dann als Anfangswerte für die DGL des gebremsten Falls.
Mit der "Event Location" kannst du "t_open" ermittlen, die Integration abbrechen und mit der neuen DGL weitermachen. Im Prinzip das gleiche wie das Ball-Beispiel auf der Matlab-Homepage.
Also ich hab mir das Beispiel aus der Matlab Hilfe zu Herzen genommen und auf meines übertragen, habe leider noch einige Fehlermeldungen! Vielleicht könntet ihr mir da noch weiterhelfen:
Code:
??? In an assignment A(I) = B, the number of elements in B and
I must be the same.
Error in ==> fallschirmdgl at 6
dy(2) = g-(rho*cw*Am)/(2*m).*y(2).^2;
Error in ==> odearguments at 111
f0 = feval(ode,t0,y0,args{:}); % ODE15I sets args{1} to yp0.
Error in ==> ode23 at 172 [neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...
Error in ==> fallschirm at 34 [t,y,te,ye,ie] = ode23(@fallschirmdgl,[tstart tfinal],y0,options);
H0 = 3000; % height of plane meter
H1 = 1500; % opening height parachut meter
Am = 0.5; % area man meter^2
Ap = 30; % area parachut meter^2
A = Am+Ap; % whole area
m = 85; % weight parachutist in kg
cw = 1.3; % cw-value
rho = 1.2; % air density kg/m^3
g = 9.81; % acceleration gravity
a1 = (cw*Am*rho)/(2*m); % air resitance before p. opens
a2 = (cw*A*rho)/(2*m); % air resitance aftzer p. opens
function[value,isterminal,direction] = events(t,y,h)
value = y(1)-h; %wenn y die Höhe zum Öffnen erreicht, dann Nulldurchgang
isterminal = 1;
direction = 0;
global g rho cw m Am;
H0 = 3000; % height of plane meter
H1 = 1500; % opening height parachut meter
Am = 0.5; % area man meter^2
Ap = 30; % area parachut meter^2
A = Am+Ap; % whole area
m = 85; % weight parachutist in kg
cw = 0;
cw_2 = 1.3; % cw-value
rho = 1.2; % air density kg/m^3
g = 9.81; % acceleration gravity
a1 = (cw*Am*rho)/(2*m); % air resitance before p. opens
a2 = (cw*A*rho)/(2*m); % air resitance aftzer p. opens
Die Variablen, die in der DGL-Funktion global sind, müssen auch in der Main-Funktion global sein ! Die Schleife, wenn man denn eine benutzt, muss natürlich nur bis 2 gehen, denn man öffnet ja nur einmal den Schirm.
Dazwischen musst du den cw-Wert auch noch ändern (geht sicher noch eleganter).
die Lösung ist noch nicht ganz richtig (pn von "robben" an mich).
Zum einen muss der Startvektor y0=[0,0] lauten, dann bekommt man für den freien Fall auch eine quadratische Abhängigkeitm, wie erwartet.
Außerdem sollte die Integration abbrechen, wenn der Springer am Boden ist; mann sollte der Event-Funktion also auch die Gesamthöhe übergeben.
Damit lautet die Main-Function:
Code:
global g rho cw m Am;
H0 = 3000; % height of plane meter
H1 = 1500; % opening height parachut meter
H=[H1,H0];
Am = 0.5; % area man meter^2
Ap = 30; % area parachut meter^2
A = Am+Ap; % whole area
m = 85; % weight parachutist in kg
cw = 0;
cw_2 = 1.3; % cw-value
rho = 1.2; % air density kg/m^3
g = 9.81; % acceleration gravity
a1 = (cw*Am*rho)/(2*m); % air resitance before p. opens
a2 = (cw*A*rho)/(2*m); % air resitance aftzer p. opens
%Phase 1 des Sprungs
tstart=0
tfinal=1000
y0=[0;0]
refine=4
%set(gca,'xlim',[030],'ylim',[0 H0]);
%box on
%hold on;
options = odeset('Events',@(t,y)events(t,y,H(i)),'OutputFcn', @odeplot,...
'OutputSel',[1,2],'Refine',refine);
% Gleichung bis interrupt lösen [t,y,te,ye,ie] = ode23(@fallschirmdgl,[tstart tfinal],y0,options);
if ~ishold
hold on
end
xlabel('time');
ylabel('height and velocity');
plot(teout,yeout(:,1),'ro',teout,yeout(:,2),'ro') title('Fallschirmsprung');
hold off
odeplot([],[],'done');
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.
RSS