作者joejoe9193 (joejoee9193)
看板MATLAB
標題[問題]function內使用尚未解出的值
時間Tue Sep 8 17:37:31 2015
一、程式說明
(1)主程式(test.m)
邊界值為Y(1/1.5)=alpha=0 Y(1)=beta=0
用shooting method去解二階ode的邊界值問題，
解ode使用的指令為ode45
(2)Function (funtest1.m)
解4條first-order initial value problems
a的值是要從判斷解出來的值運算後(也就是解出來的yv去運算)是否有大於1來設定
H=0.25;
m=1.2;
YY=yv(1,1);
YY1=yv(2,1);
si=((YY/x)^2-YY*YY1/x+(YY1)^2)^0.5
if si>1
aa=(si.^m-1)/(H*si)
elseif si<=1
aa=0
end
end
二、程式碼
test.m
alpha=0;
beta=0;
aa=1/1.5;
b=1;
% First shooting:
epsilon=1e-10;
s(1)=(beta-alpha)/(b-aa);
[X,Y]=ode45('fun1test',[aa b],[alpha;s(1);0;1]);
% iterations:
flag=0;
n=length(Y(:,1));
diff=abs(Y(n,1)-beta);
if diff<epsilon
flag=1;
end
k=1;
while flag==0,
s(k+1)=s(k)+(beta-Y(n,1))/Y(n,3);
clear X;
clear Y;
[X,Y]=ode45('fun1test',[aa b],[alpha;s(k+1);0;1]);
n=length(Y(:,1));
diff=abs(Y(n,1)-beta);
if diff<epsilon
flag=1;
end
k=k+1;
end
fun1test.m
function yv=fun1test(x,y)
global C
%C為常數
yv(1,1)=y(2,1);
yv(2,1)=1/(x*(1+a)).*C+1/(x.^2)*y(1,1)-1/x.*y(2,1);
yv(3,1)=y(4,1);
yv(4,1)=1/(x.^2)*y(3,1)-1/x*y(4,1);
三、問題
請問function裡面的a要怎麼設定?
謝謝
--
※ 發信站: 批踢踢實業坊(ptt.cc), 來自: 140.113.130.173
※ 文章網址: https://webptt.com/m.aspx?n=bbs/MATLAB/M.1441705054.A.62C.html