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

Robertson DAE Problem
 

Tommy7571

Gast

Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 18.10.2018, 15:39     Titel: Robertson DAE Problem
  Antworten mit Zitat      
Hello,

I have a problem using ode15s for solving the Robertson DAE problem.
I found the Robertson DAE problem online for matlab,

Code:

 1;
# y1' = -0.04'y1+10^4*y2*y3;
# y2' = 0.04*y1-10^4*y2*y3-3*10^7*y2^2
# y3' = 3E7*y2^2
#the conservation law is
#  y'(1) +   y'(2)   +   y'(3) = 0      
 #with initial conditions the conservation law is
#  y(1) +   y(2) + y(3) = 1        

function dy = f(t,y)
 # option 1
#  dy = [0; 0; 0];
 # dy(1)  = -0.04*y(1)+(10^4)*y(2)*y(3);                    
 # dy(2)  = 0.04*y(1)-(10^4)*y(2)*y(3)-(3*10^7)*y(2)^2;  
 # dy(3)  = 3E7*y(2)^2;
# option 2
dy = [-0.04*y(1) + 1e4*y(2).*y(3) ;  0.04*y(1) - 1e4*y(2).*y(3) - 3e7*y(2).^2 ;  y(1) + y(2) + y(3) - 1 ];    
endfunction
#t = [0, 4*logspace(-6,6)]; # ursprünglich im Code, geht auch
t = logspace(-6,6,50);

y0 = [1; 0; 0];     #initial conditions als Spaltenvektor
M = [1 0 0; 0 1 0; 0 0 0];    
options = odeset('Mass',M,'RelTol',1e-4,'AbsTol',1e-6 );
[t,y] = ode15s(@f,t,y0,options);  
y(:,2) = 1e4*y(:,2);
semilogx(t,y);
ylabel('1e4 * y(:,2)');
title('Robertson DAE problem with a Conservation Law, solved by ODE15S');
 


If I use as a result in the function the given vector (option 2) then it works but if I use option 1 - which should be just another way to express the same thing - it gives an error

"
[IDA ERROR] IDASolve
* At t = 1e-006 and h = 7.22418e-016, the corrector convergence failed repeatedly or with |h| = hmin.
error: IDASolve failed
error: called from ode15s at line 316 column 17
diff_ode15s_Robertson_DAE_problem_with_a_Conservation_Law at line 27 column 6"
(the line where the ode15 s function is used: ode15s(@f,t,y0,options)Wink
"
I would like to understand this error. Can anyone help me??

Thanks.

Regards

Tommy7571


Tommy7571

Gast

Beiträge: ---
Anmeldedatum: ---
Wohnort: ---
Version: ---
     Beitrag Verfasst am: 19.10.2018, 16:14     Titel: Robertson DAE Problem
  Antworten mit Zitat      
Hello,

I found how to run the program. The following code works:

Code:

1;
# y1' = -0.04'y1+10^4*y2*y3;
# y2' = 0.04*y1-10^4*y2*y3-3*10^7*y2^2
# y3' = 3E7*y2^2
#the conservation law is
#  y'(1) +   y'(2)   +   y'(3) = 0      
 #with initial conditions the conservation law is
#  y(1) +   y(2) + y(3) = 1        

function dy = f(t,y)
  dy = zeros(3,1); % a column vecto   #dy = [0 ; 0 ;  0];
  dy(1)=-0.04*y(1)+(1e4)*y(2).*y(3);                    
  dy(2)=0.04*y(1)-(1e4)*y(2).*y(3)-(3*10^7)*y(2).^2;   # mit Komponenten geht es auch!!!!
  dy(3)=3e7*y(2).^2;
# vorige Fehler:
# kein .*  ;  #  options = odeset('Mass',M,'RelTol',1e-4,'AbsTol',1e-6 );  absTol   falsche Länge
 
# dy = [-0.04*y(1)+1e4*y(2).*y(3);          # zuerst ging es nur damit !!
  # 0.04*y(1)-1e4*y(2).*y(3)-3e7*y(2).^2;
  # y(1)+y(2)+y(3)-1];    
endfunction
t = logspace(-6,6,50);#t = [0, 4*logspace(-6,6)]; # ursprünglich im Code, geht auch

y0 =[1; 0; 0];     #initial conditions als Spaltenvektor
M = [1 0 0; 0 1 0; 0 0 0];  
absT = [1e-6 1e-10 1e-6];  
# options = odeset('RelTol',1E-4,'AbsTol',absT); # für ode45
#options = odeset('Mass',M,'RelTol',1e-4);# ode15s  # geht nicht
 options = odeset('RelTol',1e-4,'AbsTol',absT); # für ode15s  mit Mass M geht es nicht !!
 
# [t,y] = ode45(@(t,y) f(t,y),t,y0,options); # geht nicht,  #y = lsode("f",y0,t); geht nicht
# [t,y] = ode15s(@(t,y) f(t,y),t,y0,options); # geht  
[t,y] = ode15s(@f,t,y0,options);  # ode15s(@f,t,y0,options) geht,  

y(:,2) = 1e4*y(:,2);
semilogx(t,y);
ylabel('1e4 * y(:,2)');
title('Robertson DAE problem with a Conservation Law, solved by ODE15S');
 

The result corresponds exactly to what should result.
However, I am not 100 % satisfied with this solutisince it does only work if I do not try to use the mass option in this case, which ist just the most important point.
If I replace the zero of M(3,3) by a 1, that gives the same sresult and it works!
Is there anyone knowing why the zero in the mass condition in this code leads to "
[IDA ERROR] IDASolve
At t = 1e-006 and h = 7.22418e-016, the corrector convergence failed repeatedly or with |h| = hmin.

error: IDASolve failed"

It is quite a famous example of matlab. It seems not to be an error of the initial conditions....


Thanks in advance

Thommy7571
 
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.