作者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/cn.aspx?n=bbs/MATLAB/M.1441705054.A.62C.html