Main Content

对具有两个解的 BVP 求解

此示例使用 bvp4c 和两个不同的初始估计值来求 BVP 问题的两个解。

假设有以下微分方程:

y+ey=0.

此方程具有如下边界条件:

y(0)=y(1)=0.

要在 MATLAB® 中对该方程求解,您需要先编写方程和边界条件的代码,然后为解生成合适的初始估计值,再调用边界值问题求解器 bvp4c。您可以将所需的函数作为局部函数包含在文件末尾(如本处所示),或者将它们作为单独的命名文件保存在 MATLAB 路径上的目录中。

编写方程代码

创建一个函数以编写方程代码。此函数应具有签名 dydx = bvpfun(x,y)dydx = bvpfun(x,y,parameters),其中:

  • x 是自变量。

  • y 是解(因变量)。

  • parameters 是未知参数值的向量(可选)。

求解器会自动将这些输入传递给该函数,但是变量名称决定如何编写方程代码。在本例中,可以将二阶方程重写为一阶方程组

y1=y2,

y2=-ey1.

用于编写这些方程代码的函数为

function dydx = bvpfun(x,y)
dydx = [y(2)
        -exp(y(1))];
end

编写边界条件代码

对于像此问题中的两点边界值条件,边界条件函数应该具有签名 res = bcfun(ya,yb)res = bcfun(ya,yb,parameters),具体取决于是否涉及未知参数。yayb 是求解器自动传递给函数的列向量,bcfun 返回边界条件中的残差。

对于边界条件 y(0)=y(1)=0bcfun 函数指定两个边界上的残差值都为零。在您的初始估计值中,这些残差值会强制应用于您指定给 bvpinit 的第一个和最后一个网格点。此问题的初始网格应该有 x(1) = 0x(end) = 1

function res = bcfun(ya,yb)
res = [ya(1)
       yb(1)];
end

获取初始估计值

调用 bvpinit 以生成解的初始估计值。x 的网格不需要有很多点,但是第一个点必须为 0,而最后一个点必须为 1,以正确指定边界条件。对 y 使用初始估计值,其中第一个分量为稍大于零的正数,第二个分量为零。

xmesh = linspace(0,1,5);
solinit = bvpinit(xmesh, [0.1 0]);

求解方程

使用 bvp4c 求解器求解 BVP。

sol1 = bvp4c(@bvpfun, @bcfun, solinit);

使用不同的初始估计值

使用解的不同初始估计值第二次求解 BVP。

solinit = bvpinit(xmesh, [3 0]);
sol2 = bvp4c(@bvpfun, @bcfun, solinit);

对解进行比较

绘制 bvp4c 针对不同的初始条件所计算的解。这两个解都满足规定的边界条件,但它们之间有不同行为。由于解并不始终唯一,不同行为展现出为解提供良好初始估计值的重要性。

plot(sol1.x,sol1.y(1,:),'-o',sol2.x,sol2.y(1,:),'-o')
title('BVP with Different Solutions That Depend on the Initial Guess')
xlabel('x')
ylabel('y')
legend('Solution 1','Solution 2')

Figure contains an axes object. The axes object with title BVP with Different Solutions That Depend on the Initial Guess, xlabel x, ylabel y contains 2 objects of type line. These objects represent Solution 1, Solution 2.

局部函数

此处列出了 BVP 求解器 bvp4c 为计算解而调用的局部辅助函数。您也可以将这些函数作为它们自己的文件保存在 MATLAB 路径上的目录中。

function dydx = bvpfun(x,y) % equation being solved
dydx = [y(2)
        -exp(y(1))];
end
%-------------------------------------------
function res = bcfun(ya,yb) % boundary conditions 
res = [ya(1)
       yb(1)];
end
%-------------------------------------------

另请参阅

| |

相关主题