求解具有未知参数的 BVP

打开实时脚本

以下示例说明如何使用 bvp4c 求解具有未知参数的边界值问题。

马蒂厄方程在区间 $[0, π]$ 上定义为

${y^{'}}^{'} + (λ - 2 q \cos (2 x)) y = 0$ .

当参数 $q = 5$ 时，边界条件为

$y^{'} (0) = 0$ ,

$y^{'} (π) = 0$ .

但这最多只能将 $y (x)$ 确定为一个数乘，因此需要第三个条件来指定特定解，

$y (0) = 1$ .

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

编写方程代码

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

x 是自变量。
y 是因变量。
lambda 是表示特征值的未知参数。

您可以用代换法 $y_{1} = y$ 和 $y_{2} = y^{'}$ 将马蒂厄方程写成一阶方程组，

${y_{1}}^{'} = y_{2}$ ,

${y_{2}}^{'} = - (λ - 2 q \cos (2 x)) y_{1}$ .

则对应的函数是

function dydx = mat4ode(x,y,lambda) % equation being solved
dydx = [y(2)
      -(lambda - 2*q*cos(2*x))*y(1)];
end

注意：所有函数都作为局部函数包含在示例的末尾。

编写边界条件代码

现在，编写一个函数，该函数返回在边界点处的边界条件的残差值。此函数应具有签名 res = mat4bc(ya,yb,lambda)，其中：

ya 是在区间 $[a, b]$ 开始处的边界条件的值。
yb 是在区间 $[a, b]$ 结束处的边界条件的值。
lambda 是表示特征值的未知参数。

此问题在区间 $[0, π]$ 内有三个边界条件。要计算残差值，您需要将边界条件设置为 $g (x, y) = 0$ 形式。在此形式中，边界条件是

$y^{'} (0) = 0$ ,

$y^{'} (π) = 0$ ,

$y (0) - 1 = 0$ .

则对应的函数是

function res = mat4bc(ya,yb,lambda) % boundary conditions
res = [ya(2)
       yb(2)
       ya(1)-1];
end

创建初始估计值

最后，创建解的初始估计值。您必须对两个解分量 $y_{1} = y (x)$ 和 $y_{2} = y^{'} (x)$ 以及未知参数 $λ$ 提供初始估计值。只有接近初始估计值的特征值和特征函数才可以被计算。为了增大计算的特征函数对应于第四个特征值的可能性，您应选择具有正确的定量行为的初始估计值。

对于此问题，余弦函数满足三个边界条件，因此有助于提供较好的初始估计值。使用返回 $y_{1}$ 和 $y_{2}$ 的估计值的函数，编写 $y$ 的初始估计值的代码。

function yinit = mat4init(x) % initial guess function
yinit = [cos(4*x)
        -4*sin(4*x)];
end

使用区间为 $[0, π]$ 的 10 点网格、初始估计值函数以及 $λ$ 的估计值 15 调用 bvpinit。

lambda = 15;
solinit = bvpinit(linspace(0,pi,10),@mat4init,lambda);

求解方程

使用 ODE 函数、边界条件函数和初始估计值调用 bvp4c。

sol = bvp4c(@mat4ode, @mat4bc, solinit);

参数值

打印 bvp4c 求得的未知参数 $λ$ 的值。此值是马蒂厄方程的第四个特征值 ( $q = 5$ )。

fprintf('Fourth eigenvalue is approximately %7.3f.\n',...
            sol.parameters)

Fourth eigenvalue is approximately  17.097.

对解进行绘图

使用 deval 计算 bvp4c 在区间 $[0, π]$ 中的 100 个点处计算的解。

xint = linspace(0,pi);
Sxint = deval(sol,xint);

对两个解分量进行绘图。绘图显示了与第四个特征值 $λ_{4} = 17.097$ 相关联的特征函数（及其导数）。

plot(xint,Sxint)
axis([0 pi -4 4])
title('Eigenfunction of Mathieu''s Equation.') 
xlabel('x')
ylabel('y')
legend('y','y''')

Figure contains an axes object. The axes object with title Eigenfunction of Mathieu's Equation., xlabel x, ylabel y contains 2 objects of type line. These objects represent y, y'.

局部函数

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

function dydx = mat4ode(x,y,lambda) % equation being solved
q = 5;
dydx = [y(2)
      -(lambda - 2*q*cos(2*x))*y(1)];
end
%-------------------------------------------
function res = mat4bc(ya,yb,lambda) % boundary conditions
res = [ya(2)
       yb(2)
       ya(1)-1];
end
%-------------------------------------------
function yinit = mat4init(x) % initial guess function
yinit = [cos(4*x)
        -4*sin(4*x)];
end
%-------------------------------------------

另请参阅

bvp4c | bvp5c | bvpinit