[Matlab科学计算] 开发解二阶常系数常微分非齐次方程u′′ + pu' + qu = f(x)源代码，a<=x<=b

以下是一个解二阶常系数常微分非齐次方程的 MATLAB 源代码：

function [x, y] = solve_ode(p, q, f, a, b, y0, y1, n)
% p, q, f: 常微分方程中的常数函数
% a, b: 自变量的取值范围
% y0, y1: 初值条件
% n: 等分数

h = (b - a) / (n - 1);
A = zeros(n, n);
b = zeros(n, 1);

% 构造系数矩阵和右端向量
for i = 2 : n - 1
    xi = a + (i - 1) * h;
    A(i, i - 1) = -1;
    A(i, i) = 2 + h^2 * q(xi);
    A(i, i + 1) = -1;
    b(i) = h^2 * f(xi);
end

% 处理边界
A(1, 1) = 1;
A(n, n) = 1;
b(1) = y0;
b(n) = y1;

% 解线性方程组
y = A \ b;

% 返回自变量和函数值
x = linspace(a, b, n)';

使用示例：

% 解 u'' + 2u' + u = cos(x), 0 <= x <= pi, u(0) = 0, u(pi) = 1
p = @(x) 2;
q = @(x) 1;
f = @(x) cos(x);
a = 0;
b = pi;
y0 = 0;
y1 = 1;
n = 100;

[x, y] = solve_ode(p, q, f, a, b, y0, y1, n);
plot(x, y);

这个例子解决了方程 $u'' + 2u' + u = \cos(x)$ 在 $0 \leq x \leq \pi$ 上的解，且边界条件为 $u(0) = 0$ 和 $u(\pi) = 1$。