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

Anfang: Funktion plotten und DGL System lösen

 

sutoppu
Forum-Anfänger

Forum-Anfänger


Beiträge: 13
Anmeldedatum: 26.04.11
Wohnort: KÖLN
Version: ---
     Beitrag Verfasst am: 26.04.2011, 07:16     Titel: Anfang: Funktion plotten und DGL System lösen
  Antworten mit Zitat      
Guten morgen Smile

ich bin grad ein wenig am Verzweifeln... Ich muss eine Funktion f plotten und danach ein DGL System lösen (die Funktion ist teil der DGLn). Ich komme allerdings gar nicht zurecht (mein erster Versuch alleine irgendwas mit Matlab zu machen).

mein Code sieht so aus (erstmal nur die Funkton):


Code:

t_p = 0.3;
t_r = 0.175;
t_c = 0.15;
alpha = 2;
t_b = t_p.*(1-((t_r./t_c).^(alpha./(alpha-1))).*((exp(-(t_p./t_r).^(alpha)))./(1-exp(t_p./t_c).^(alpha))).^(1/(alpha-1)));
h = 0.0001;
t_end=1;
N = t_end/h;
g3 = zeros(N+1,1);
g4 = zeros(N+1,1);
f = zeros(N+2,1);
t = zeros(N+1,1);
i=1;

while t(i) < t_end;
    g3(i+1) = ((1-exp(-(t(i)./t_c).^(alpha)))*(exp(-((t(i)-t_b)/t_r).^(alpha))))./((1-exp(-(t_p./t_c).^(alpha)))*(exp(-((t_p-t_b)/t_r).^(alpha))));
    g4(i+1)=(1-exp(-(t(i)./t_c).^(alpha)))./((1-exp(-(t_p./t_c).^(alpha)))*(exp(-((t_p-t_b)/t_r).^(alpha))));
      if t(i)<= 0.3;
        f(i)= g4(i);
      else
          f(i)=g3(i);
     
       end;
t(i+1) = t(i)+h;
    i=i+1;
end;


ist wahrscheinlich sehr umständlich... ich wusste nicht, ob man eine Funktion aus zwei Teilen mit einem normalen Befehl plotte kann... in der if-abfrage sollte eigentlich t_b stehen, aber dann bekomm ich eine Fehlermeldung, also habe ich den Wert für t_b ausgerechnet... naja Smile Das Bildchen sieht der Lösung aus dem Buch auch ziemlich ähnlich!

Kann mir jemand sagen, ob das so ok ist?

Funktion.jpg
 Beschreibung:
Funktion f

Download
 Dateiname:  Funktion.jpg
 Dateigröße:  49.65 KB
 Heruntergeladen:  732 mal
Private Nachricht senden Benutzer-Profile anzeigen


sutoppu
Themenstarter

Forum-Anfänger

Forum-Anfänger


Beiträge: 13
Anmeldedatum: 26.04.11
Wohnort: KÖLN
Version: ---
     Beitrag Verfasst am: 26.04.2011, 08:51     Titel: weiter gehts
  Antworten mit Zitat      
Wenn ich die Funktion f habe, dann muss ich ein System lösen. Ich wollte es erstmal mit nem expl. Euler Versuchen... Das System besteht aus zwei Dgln.
Ich hab versucht als erstes das System zu lösen, aber das klappt irgendwie schon nicht...

Code:

t_p = 0.3;
t_r = 0.175;
t_c = 0.15;
alpha = 2;
t_b = t_p.*(1-((t_r./t_c).^(alpha./(alpha-1))).*((exp(-(t_p./t_r).^(alpha)))./(1-exp(t_p./t_c).^(alpha))).^(1/(alpha-1)));
R0=0.08;
h = 0.0001;
t_end=1;
N = t_end/h;
pr=10;
Rin=0.001;
Cs=2.75;
Rs=1.0;
a=0.0007;
b=5;
c=1.6;
d=1;

t = zeros(N+1,1);
x1 = zeros(N+1,1);
x2 = zeros(N+1,1);
x1(1) = 72;
x2(1) = 125;
t(1) = 0;
i = 1;
while t(i) < t_end
    f1=(-(((Rs+R0)/(R0*Rs*Cs)))*x1(i))+(1/(R0*Cs))*(a*((x2(i)-b)^2)+((c*x2(i)-d)*f(i)));
    f2=(1/R0)*x1(i)-((1/R0)*((a*(x2(i)-b).^2)+((c*x2(i)-d)).*f(i)));
    x1(i+1) = f1*h + x1(i);
    x2(i+1) = f2*h + x2(i);
    t(i+1) = t(i) + h;
    i=i+1;
end
pv=[(a*((x2-b).^2))+((c*x2-d).*f)];
p_iso = a*(125-b).^2+(c*125-d)*f;
pa=x1-R0*f2;
plot(t,pv,t,p_iso,t,pa);
 


p_iso stimmt so. Aber pv und pa sind falsch und ich sehe nicht warum... Meine Lösung von V (im Code f2, bzw x2) und ps (im Code f1 bzw x2) sind glaube ich schon nicht richtig.
Und pa kann man so glaube ich auch nicht berechnen, aber mir ist bisher nix besseres eingefallen....

Ich hab das Bild, das rauskommen soll mal angehängt... vllt liegts aber auch schon an der Funktion f aus dem ersten Teil?

Ich hoffe jemand hat eine Idee Smile

Lösung.jpg
 Beschreibung:

Download
 Dateiname:  Lösung.jpg
 Dateigröße:  98.89 KB
 Heruntergeladen:  731 mal
system.jpg
 Beschreibung:

Download
 Dateiname:  system.jpg
 Dateigröße:  36.18 KB
 Heruntergeladen:  722 mal
Private Nachricht senden Benutzer-Profile anzeigen
 
Thomas84
Forum-Meister

Forum-Meister


Beiträge: 546
Anmeldedatum: 10.02.10
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 27.04.2011, 05:53     Titel:
  Antworten mit Zitat      
Ich würde vorschlagen du schaust dir mal die DGL-Löser von Matlab an (ode45). Es ist glaub ich noch einiges verkehrt bei deinem code. Du benötigst z.B. pv innerhalb der Schleife berechnest es aber erst nachher. Außerdem ist

Code:

x2(i+1) = f2*h + x2(i);
 


sicherlich nicht richtig!?
Private Nachricht senden Benutzer-Profile anzeigen
 
sutoppu
Themenstarter

Forum-Anfänger

Forum-Anfänger


Beiträge: 13
Anmeldedatum: 26.04.11
Wohnort: KÖLN
Version: ---
     Beitrag Verfasst am: 27.04.2011, 17:16     Titel:
  Antworten mit Zitat      
Ich habe im Code noch Zeichenfehler gefunden und das Problem mit pa gelöst Smile ich musste dV/dt > 0 ausgeben. Damit haben sich natürlich die anderen Lösungen stark verändert. Und meine Lösung sieht schon ganz gut aus. Aber ich meine dass f(i) noch nicht ganz richtig ist, für alpha <2 sieht die Lösung nicht aus wie eine Parabel...
Pv habe ich einmal in der Schleife als auch danach berechnet, aber die Lösung ist dieselbe.
Dafür musste ich aber pa innerhalb der Schleife bestimmen

ich werd mir ode45 mal angucken! Aber warum sollte die Vorschrift nicht stimmen? Das sind eigentlich die einzigen zwei Zeilen bei denen ich mir sicher war Wink
Private Nachricht senden Benutzer-Profile anzeigen
 
Thomas84
Forum-Meister

Forum-Meister


Beiträge: 546
Anmeldedatum: 10.02.10
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 28.04.2011, 12:01     Titel:
  Antworten mit Zitat      
Ich glaub ich hatte den Code noch nicht ganz durchschaut. Die Verwendung von ode45 und ein paar Kommentare würden denn Code sicherlich übersichtlicher machen.
Private Nachricht senden Benutzer-Profile anzeigen
 
sutoppu
Themenstarter

Forum-Anfänger

Forum-Anfänger


Beiträge: 13
Anmeldedatum: 26.04.11
Wohnort: KÖLN
Version: ---
     Beitrag Verfasst am: 01.05.2011, 12:10     Titel: neues Problem
  Antworten mit Zitat      
Hallo Smile

ich hab das System gelöst! Danke, dass ihr euch mein Problem angesehen habt Smile das nächste mal, werde ich den Code auch ordentlicher aufschreiben!

Jetzt habe ich aber ein neues Problem:

In einem weiteren System habe ich eine Zeitverzögerung:

Code:

p_s= (f(i)-k1*(-v_w(i)) + k2*((-v_w(i-kappa*i)).^2))));
 


mit kappa = 0.45

wenn ich das so implementiere, dann bekomme ich einen Fehler, dass man nur integer eingeben darf... wie kann man das beheben?
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: 01.05.2011, 12:46     Titel:
  Antworten mit Zitat      
Hallo,

so etwas wird als "Delay Differential Equation" bezeichnet.
Fertige Löser sind dde23 oder ddesd, siehe
http://www.mathworks.com/access/hel.....chdoc/math/f1-663833.html

Wenn du das selber (mit fester Schrittweite) löst, dann solltest du die Simulationsschrittweite so wählen, dass die Verzögerung ein vielfaches ist, z.B. hier dt = 0.15. Dann entspricht kappa = 0.45 einer Verzögerung um 3 Zeitschritte.

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

Forum-Anfänger

Forum-Anfänger


Beiträge: 13
Anmeldedatum: 26.04.11
Wohnort: KÖLN
Version: ---
     Beitrag Verfasst am: 01.05.2011, 12:56     Titel:
  Antworten mit Zitat      
Wow, danke für die schnelle antwort Smile

Ich werds direkt mal probieren Smile
Private Nachricht senden Benutzer-Profile anzeigen
 
sutoppu
Themenstarter

Forum-Anfänger

Forum-Anfänger


Beiträge: 13
Anmeldedatum: 26.04.11
Wohnort: KÖLN
Version: ---
     Beitrag Verfasst am: 01.05.2011, 15:42     Titel:
  Antworten mit Zitat      
Hallo!

Ich habs jetzt erstmal mit fester Schrittweite probiert, es hat soweit auch funktioniert. Problem ist jetzt nur, dass meine Anfangswerte nicht mehr passen Sad also geht jetzt die Fehlersuche los...

Ich probiers erstmal alleine, sonst meld ich mich nochmal Smile

DANKE!
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.