用matlab求y=f(x)，y''=y^2的数值解，边界问题

Question

待解方程为：y=f(x), y''=y^2，即y的二阶导等于y的平方
边界条件：y(0)=1, y'(100)=0
用matlab函数ode45求解
function dy=test(x,y)
dy=zeros(2,1);
dy(1)=y(2);
dy(2)=[y(1)]^1;
end
[X,Y]=ode45('test',[0 100],[1 ?]); 这个地方就不会写了
看教程这个‘？’地方应该是填y'(0)的值，由于不知道y'(0)的值，只知道y'(100)=0
求大神指点迷津，谢谢。

Answer 1

分析

这是个边值问题（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''');

 

结果

Answer 2

可以提供两种方法：
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