文章目录
一、思路
霍纳方法(或秦九韶方法)是一种高效的计算
n
n
n次多项式
P
n
(
u
)
P_n(u)
Pn(u)值的算法,其具体推导见博文:多项式求值。由于对n次多项式求d阶导数,还是多项式,因而,很自然地,n次多项式的任意阶导数计算,也可以采用霍纳方法。关键在于如何确定d阶多项式的系数?
不失一般性,这里以3次多项式为例,观察其d阶导数的多项式系数的变化规律:
(1)3次多项式1阶导数的多项式系数相对于0阶导数的多项式系数,需要更新最后2个。最后一个系数乘以次数3,得到
3
a
3
3a_3
3a3,倒数第二个系数乘以2(即次数减1),得到
2
a
2
2a_2
2a2。
(2)3次多项式2阶导数的多项式系数相对于1阶导数的多项式系数,需要更新最后1个。最后一个系数乘以2(即次数减1),得到
6
a
3
6a_3
6a3。
(3)3次多项式3阶导数的多项式系数相对于2阶导数的多项式系数,不需要更新。
更一般地,可以总结得到:
n次多项式d阶导数的多项式系数相对于d-1阶导数的多项式系数,需要更新最后n-d个,d阶导数的第j个系数=d-1阶导数的第j个系数 *(j - d)
二、MATLAB代码
%{
Function: calc_poly_ders
Description: 霍纳方法(或秦九韶方法)计算n次多项式直到d阶导数
Input: 多项式结构体,自变量u,导数阶数d
Output: n次多项式0到d阶导数ders
Author: Marc Pony(marc_pony@163.com)
%}
function ders = calc_poly_ders(poly, u, d)
%赋值给临时变量,方便书写
a = poly.a;
us = poly.us;
order = poly.order;
ders = zeros(1, d + 1);
b = a; %赋值给临时变量
%计算第0阶导数
ders(1) = b(order + 1);
for i = order : -1 : 1
ders(1) = (u - us) * ders(1) + b(i);
end
%计算第1到d阶导数
minD = min([d, order]);
for k = 1 : minD
for j = order + 1 : -1 : k + 2
b(j) = (j - k) * b(j);
end
ders(k + 1) = b(order + 1);
for i = order : -1 : k + 1
ders(k + 1) = (u - us) * ders(k + 1) + b(i);
end
end
end
clc
clear
poly = struct();
poly.order = 5;
poly.us = rand();
poly.a = rand(1, poly.order + 1);
u = rand();
d = 3;
pos=poly.a(1) + poly.a(2)*(u - poly.us) + poly.a(3)*(u - poly.us)^2 + poly.a(4)*(u - poly.us)^3 + poly.a(5)*(u - poly.us)^4 + poly.a(6)*(u - poly.us)^5;
vel=poly.a(2) + 2*poly.a(3)*(u - poly.us) + 3*poly.a(4)*(u - poly.us)^2 + 4*poly.a(5)*(u - poly.us)^3 + 5*poly.a(6)*(u - poly.us)^4;
acc=2*poly.a(3) + 6*poly.a(4)*(u - poly.us) + 12*poly.a(5)*(u - poly.us)^2 + 20*poly.a(6)*(u - poly.us)^3;
jerk=6*poly.a(4) + 24*poly.a(5)*(u - poly.us) + 60*poly.a(6)*(u - poly.us)^2;
ders_ = [pos, vel, acc, jerk]
ders = calc_poly_ders(poly, u, d)