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

ODE45 mit Events
 

manuel-690
Forum-Newbie

Forum-Newbie

Beiträge: 3
Anmeldedatum: 07.05.13
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 07.05.2013, 22:35     Titel: ODE45 mit Events
  Antworten mit Zitat      
Tag miteinander,

ich habe ein kleines Problem bezüglich der Event Funktion, die einen ODE45 Integrator überwachen soll.

Bei meinem eindimensionalen Modell handelt es sich um einen springenden Ball. Dieser wird lediglich aus einer Anfangshöhe fallen gelassen. Der Sprung am Boden soll durch eine Feder simuliert werden.
Demnach überwacht die Event Funktion den Nulldurchgang der Ballunterseite, damit die Normalkraft angeschaltet und der Ball wieder nach oben beschleunigt werden kann.

Mein Problem liegt darin, dass der Ball nach einem Bounce nie wieder die "Eintauchgeschwindigkeit" erreicht und somit bei jedem Bounce an Höhe verliert.

Hier beispielsweise ein Auszug aus der Matrix "xout" [y ; dy] mit einem Ballradius = 2 m.


1.999999999999940 -18.292876463820740 <--- Eintauchen bei y = 2m
1.819618375192020 -17.482297359361933
1.655794485106991 -14.944457864007815
1.525355011360619 -10.920471635168091
1.441353460247758 -5.798291899618427
1.411558004825795 -0.109137511925919
1.439129665944135 5.586766859952144
1.521495219963670 10.727157567665689
1.650221266115319 14.808449451642741
1.812695913775857 17.415233370706240
1.992520389371559 18.258815842964697
2.000000000000063 18.260927418190910 <--- Auftauchen bei y = 2 m

Hier noch der Code:

Code:
%
clc
clear all
close all

global g R m k forceon

format long

%% initital conditions

g = 9.81;
R = 2 ;
m = 10 ;
k = 10000;


ts = 0;
tf = 100.0;

x0=[20;0];       %x= = [y ; dy]
tout = ts;
xout = x0.';
teout = [];
xeout = [];
ieout = 0;
FO = [0];
forceon = 0;

%% integrator

opt = odeset('Events',@event_1D_event,'AbsTol',1e-6,'RelTol',1e-3);

for i = 1:50
    if ts < (tf-0.01)
       
        t0 = ts:0.01:tf;
        [t, x,te,xe,ie]=ode45(@DGL_1D_event, t0 ,x0,opt);
        if ~ishold
            hold on
        end
        nt = length(t);
        tout = [tout; t(2:nt)];
        xout = [xout; x(2:nt,:)];
        teout = [teout;te];
        xeout = [xeout;xe];
        ieout = [ieout;ie];
       
       
        for j = 2:nt
        FO = [FO;forceon];
        end
       
         x0 = CallBack_1D_event(nt,x);
       
    ts = t(nt);
    end
end;


%% plot

figure(1)
plot(tout,xout(:,1),tout,FO);
grid on
legend('y')
xlabel('t [s]')
ylabel('y [m]')

----------------------------------------------------------------------

function [value,isterminal,direction] = event_1D_event(t,x)

global R

value = (x(1)-R);
isterminal = 1;
direction = 0;

end

----------------------------------------------------------------------

function dx=DGL_1D_event(t,x)

global g R m k forceon
     
   
%% DGL
 
if forceon == 1
    beschly = ((k*(abs(x(1)-R)))/m)-g;
else
    beschly = - g;
end

dx=[x(2);beschly];


----------------------------------------------------------------------

function newx = CallBack_1D_event(nt,x)

global g R m k forceon

   newx = [x(nt,1);x(nt,2)];
   
if (x(nt,1)-R) <= 0 && (forceon == 0)
    forceon = 1;
    return
end

if ((x(nt,1)-R) >= 0 ) && forceon == 1
    forceon = 0;
    return
end
       
end
 


Ich hoffe das war soweit verständlich. Vielen Dank im Voraus.
Private Nachricht senden Benutzer-Profile anzeigen


Harald
Forum-Meister

Forum-Meister

Beiträge: 24.501
Anmeldedatum: 26.03.09
Wohnort: Nähe München
Version: ab 2017b
     Beitrag Verfasst am: 07.05.2013, 23:09     Titel:
  Antworten mit Zitat      
Hallo,

wenn ich das richtig sehe, musst du nur die Toleranzen enger machen.
Code:
opt = odeset('Events',@event_1D_event,'AbsTol',1e-9,'RelTol',1e-9);


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

 Forum-Newbie

Forum-Newbie

Beiträge: 3
Anmeldedatum: 07.05.13
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 07.05.2013, 23:15     Titel:
  Antworten mit Zitat      
Tatsächlich das war's schon.
Vielen Dank für die schnelle Antwort.
Ich dachte das hätte ich schon probiert und habe eine Fehlermeldung bekommen aber jetzt funktioniert es. Smile Da ist die Fehlertoleranz nur noch 10^-5.
Super!
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 - 2025 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.