分析
这是个边值问题(BVP),不是初值问题。求解边值问题不能用ode系列函数,需要用专门的求解器。下面给你参考代码,涉及到的函数主要有bvpinit、bvp4c、deval,请自行参阅相关函数的说明。
参考代码
dydx = inline('[y(2); y(1)^2]', 'x', 'y');
bc = inline('[ya(1)-1; yb(2)]', 'ya', 'yb');
sol = bvpinit(linspace(0,100,50),[0 0]);
sol = bvp4c(dydx, bc, sol);
x = linspace(0,100);
y = deval(sol,x);
subplot 211
plot(x,y(1,:));
xlabel('Time (s)'); ylabel('y');
subplot 212
plot(x,y(2,:));
xlabel('Time (s)'); ylabel('y''');
结果
可以提供两种方法:
1:迭代法,通过自变量步长推进求解,有一定的算法。
2:MATLAB符号运算的自带函数dsolve,可以求出解符号表达式,用自变量的域代替就行了。
y=dsolve('D2y-3*Dy=x^2','Dy(0)=1','y(1)=0','x');
如果要求-10到10之间的解(Y值),可令步长为0.01
x=(-10:0.01:10);
y的解:
y=subs(y,'x',x)
欧拉法就是差分法,首先要给定定义域,然后离散定义域,给定边界条件,然后可以用显示推进或隐式迭代出每一个步长的Y值。一般可以先降阶,将二阶降到一阶。显式一班有如下形式:
for i=1:n
y 1= y + h*feval(f,t,y); % evaluate slope= feval(f,t,y);
t = t + h;
tF(i+1,:) = t; % store new values in arrays tF, yF
yF(i+1,:) = y.'; % Note: y values are stored as ROW of array yF
end
隐式有
for i=1:n
x=?
while abs(temp-y)>0.001
y=f(x)
temp=y
end
end