フィルターのクリア

loop for with ode45

2 ビュー (過去 30 日間)
tony
tony 2013 年 1 月 23 日
hi , i've got a problem..
I want to solve a EDO , so i used ODE45, i've got two function :
*function dy=f_iso(t,y)
global Q b K n E sigma0 epoint
if abs(y(1))-Q*(1-exp(-b*y(2)))-sigma0<0
dy(2)=0
else
dy(2)=((abs(y(1))-Q*(1-exp(-b*y(2)))-sigma0)/K).^n
end
dy(3)=sign(y(1))*dy(2);
dy(1)=E*(epoint-dy(3));
dy=dy';
end*
And my main function :
[Tn,N]=ode45('f_iso',[0 10],[0 0 0]);
sigma1=N(:,1);
ep1=N(:,3);
e1=ep1+sigma1./E;
It's the first court , then i did the second : [Th,h]=ode45('f_iso',[10 30],[N(89,:]);
sigma2=h(:,1);
ep2=h(:,3);
e2=ep2+sigma2./E;
plot(e1,sigma1,e2,sigma2)
and i want to do it 10 times , so i did a loop for :
k=3
for i=1:k
Ti=10
Tf=30
epoint=(-1)^i*0.001;
C=N(length(N),:)
[T,b]=ode45('f_iso',[Ti Tf],[C])
Ti=Ti+20
Tf=Tf+20
C=b(length(b),:);
end
but i doest work and i don't know how to use plot with a loop for...
tranks :)
  1 件のコメント
Jan
Jan 2013 年 1 月 24 日
Please format your code properly. Currently it is unnecessarily hard to read your question.
You did not answer my questions for clarifications in your related question http://www.mathworks.com/matlabcentral/answers/59582-ode45-errer-input-argument-y-is-undefined. There I've explained already, that ODE45 is not sufficient to solve a discontinous system. I strongly suggest to consider this.
Please explain what "it does not work" exactly means: Do you get an error message or does the result differ from your expectations.

サインインしてコメントする。

回答 (1 件)

Jan
Jan 2013 年 1 月 24 日
編集済み: Jan 2013 年 1 月 24 日
You reset Ti, Tf and C in each iteration to the same values. Perhaps you want this:
global epoint % !!!
Ti = 10;
Tf = 30;
C = N(size(N, 1), :);
for i=1:k
epoint = (-1)^i*0.001; % Not used at all?
[T, b] = ode45(@f_iso, [Ti Tf], C);
Ti = Tf
Tf = Ti+20
C = b(size(b, 1), :);
end
Using global variables to define parameters for ODE45 is not smart, see http://www.mathworks.de/matlabcentral/answers/1971.

カテゴリ

Help Center および File ExchangeLoops and Conditional Statements についてさらに検索

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by