필터 지우기
필터 지우기

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 CenterFile Exchange에서 Loops and Conditional Statements에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by