14.2.1梯形法则
完整代码
trap函数
14.2.2 辛普森法则
14.2.3 quad
这其实相当简单，它依赖于一个事实，即函数 在 和 之间的定积分等于由 轴和垂直线 和 所围成的面积。因此，求解积分的所有数值方法都涉及以更多或更少巧妙的方式估计 下的面积。
梯形法则是编程实现起来相当简单的数值积分方法。它将函数 下的面积分成宽度为 的垂直面板，称为步长。如果有 个这样的面板，那么 ，即 。如果我们连接相邻面板与 的交点，我们可以估计 下的面积，其为最终得到的梯形的面积之和。如果我们将这个对积分的近似值称为 ，则
其中 。上式为梯形法则，给出了积分的估计.
这是一个用于实现梯形法则的函数：
function y = trap( fn, a, b, h )
n = (b-a)/h;
x = a + [1:n-1]*h;
y = sum(feval(fn, x));
y = h/2*(feval(fn, a) + feval(fn, b) + 2*y);
输入参数：
计算步长数量 n：根据输入的起始点和终点，以及步长 h，计算步长数量 n，即 (b-a)/h。
计算梯形面板上的 x 值：创建一个长度为 n-1 的向量 x，用于存储梯形面板上的 x 值。通过计算 x = a + [1:n-1]*h 得到每个面板的 x 值。
计算面板上函数值的总和：使用 feval 函数调用输入的函数 fn，在每个面板的 x 值处计算函数值，并将这些函数值相加。
计算积分的近似值：使用梯形法则的数值公式计算积分的近似值。首先，计算首尾两个点的函数值，然后将其与面板上函数值的总和相加，并乘以步长 h 的一半。
返回积分的近似值：将计算得到的近似值 y 返回给调用者。
需要注意的是：
.M 文件实现中使用适当的数组操作符，以便能够进行向量化计算。通过调用这个 trap 函数，可以实现对函数的数值积分，使用梯形法则进行面积估计。
计算函数 在区间 0 到 4 上的积分。在函数 .M 文件中，需要将 定义为 x.^3。可以按照以下方式调用 trap 函数来计算积分估计值：
h = 0.1
s = trap(@f, 0, 4, h)
这将使用步长 来计算积分的近似值为 64.04。如果将步长改为 ，则积分的近似值为 64.0004（精确积分值为 64）。可以观察到，随着步长的减小，积分的估计值变得更加准确。
这个例子假设 是一个连续函数，可以在任意 处进行评估。在实际应用中，函数可能是在离散点上定义的，作为实验结果的一部分提供。例如，一个物体的速度 可能每隔一段时间测量一次，人们可能希望通过速度-时间曲线下的面积来估计行驶距离。在这种情况下，trap 函数需要修改，通过替换 fn 为函数值的向量来处理离散数据。这个修改可以作为一个有趣的练习。另外，你也可以使用 MATLAB 的 interp1 函数来对数据进行插值。
% 定义函数 f(x)
function y = f(x)
y = x.^3;
end
% 梯形法则计算积分估计值的函数
function y = trap(fn, a, b, h)
n = (b - a) / h;
x = a + (1:n-1) * h;
y = sum(feval(fn, x));
y = h/2 * (feval(fn, a) + feval(fn, b) + 2*y);
end
% 主程序
h = 0.1; % 步长
a = 0; % 区间下限
b = 4; % 区间上限
% 调用 trap 函数计算积分估计值
s = trap(@f, a, b, h);
disp(s); % 打印积分估计值
trap 函数是用来计算梯形法则的积分估计值的。
function y = trap(fn, a, b, h)
fn：函数句柄（function handle），指定要进行积分的函数。函数句柄用 @ 符号后跟函数名的方式来表示。a：积分区间的下限。b：积分区间的上限。h：步长，定义每个梯形面板的宽度。函数主体：
n = (b - a) / h;
x = a + (1:n-1) * h;
y = sum(feval(fn, x));
y = h/2 * (feval(fn, a) + feval(fn, b) + 2*y);
h 计算出面板的数量 n。[a, b] 区间上生成每个面板上的 x 值，存储在向量 x 中。feval 函数，在每个面板的 x 值处计算函数值，并将这些值相加得到面板上函数值的总和。y。返回值：
y：计算出的积分的近似值。使用时，需要提供要进行积分的函数的函数句柄，设置积分区间的下限 a 和上限 b，以及所需的步长 h，并调用 trap 函数即可计算积分的近似值。
示例用法：
h = 0.1; % 步长
a = 0; % 区间下限
b = 4; % 区间上限
% 调用 trap 函数计算积分估计值
s = trap(@f, a, b, h);
以上示例中，trap 函数被用来计算函数 f(x) = x^3 在区间 [0, 4] 上的积分估计值，步长为 0.1。计算结果存储在变量 s 中。根据需求，可以对 a、b 和 h 进行调整以适应不同的积分计算。
Simpson's rule（辛普森规则）是一种数值积分方法，比梯形法则更准确，通常在尝试其他更复杂的方法之前应使用它。辛普森规则将要积分的函数 的区域分成垂直条带，但与用直线连接点 不同的是，每个连续的三个点被拟合成一个抛物线。为了确保总是有偶数个面板，通常选择步长 ，使得面板的数量为 ，即 。使用与之前相同的符号，辛普森规则通过以下公式估计积分值：
将这个公式编码为函数.M文件留做练习。
辛普森规则具有一个有趣的额外优势：它能精确地积分三次多项式（可以证明）。因此，如果你尝试在任意限制下应用辛普森规则来计算 的积分，你会发现令人惊讶的是，它给出的结果与精确的数学解相同。
下面是将辛普森规则公式编码为函数.m文件的示例代码：
function S = simpsons_rule(f, a, b, n)
h = (b - a) / (2 * n);
% 计算奇数索引和偶数索引的函数值总和
sum_odd = 0;
sum_even = 0;
for i = 1:n-1
x = a + 2*i*h; % 奇数索引处的x值：a + 2i*h
sum_odd = sum_odd + feval(f, x);
x = a + (2*i+1)*h; % 偶数索引处的x值：a + (2i+1)*h
sum_even = sum_even + feval(f, x);
end
% 计算积分估计值
S = (h/3) * (feval(f, a) + feval(f, b) + 2*sum_odd + 4*sum_even);
end
上述代码定义了一个名为 simpsons_rule 的函数。它接受以下参数：
f：积分函数的函数句柄。a：积分区间的下限。b：积分区间的上限。n：面板的数量，这决定了步长（h）的大小。函数的主体中，首先根据步长公式计算步长 h。然后使用循环计算奇数索引和偶数索引处的函数值总和。最后，使用辛普森规则的公式计算积分的近似值，并将结果返回。
使用此函数，可以进行以下示例调用：
f = @(x) x.^3; % 定义积分函数
a = 0; % 区间下限
b = 2; % 区间上限
n = 10; % 面板的数量
S = simpsons_rule(f, a, b, n); % 调用辛普森规则函数计算积分
disp(S); % 输出积分的近似值
以上示例中，我们使用积分函数 f(x) = x^3 在区间 [0, 2] 上应用辛普森规则来计算积分的近似值。我们选择 10 个面板来进行估计。计算结果存储在变量 S 中，并通过 disp 函数输出。根据需要，你可以调整函数的参数来适应不同的积分计算。
quad函数是MATLAB中用于数值积分（quadrature）的函数。它可以用于计算定积分，即对一个给定函数在指定区间上的积分值。quad函数的语法如下：
Q = quad(fun, a, b)
Q = quad(fun, a, b, tol)
Q = quad(fun, a, b, tol, args)
参数说明：
fun：代表要进行积分的函数，可以是匿名函数、函数句柄或函数名称。a 和 b：代表积分区间的上下界。tol：代表可接受的积分误差，默认为1e-6。较小的误差值会增加计算时间，但提供更精确的结果。args：可选参数，传递给 fun 函数的附加参数。quad函数使用的是自适应数值积分算法，它根据积分区间内的函数变化情况自动调整积分的步长。该算法会将积分区间分成多个小段，并在每个小段上进行数值积分，然后将结果相加得到最终的积分值。
quad函数返回的结果是一个数值，代表了对指定函数在指定区间上的数值积分值。请注意，quad函数仅适用于连续的函数。对于不连续或具有特殊性质的函数，可能需要使用其他数值积分方法或手动处理。
使用quad函数时，可以通过调整容差(tol)参数来控制积分的精确度和计算效率。较小的容差值将产生更准确的积分结果，但可能需要更长的计算时间。
联系客服
