n次多项式的任意阶导数计算(附代码)

一、思路

  霍纳方法(或秦九韶方法)是一种高效的计算$n$次多项式$P_n(u)$值的算法，其具体推导见博文：多项式求值。由于对n次多项式求d阶导数，还是多项式，因而，很自然地，n次多项式的任意阶导数计算，也可以采用霍纳方法。关键在于如何确定d阶多项式的系数？
  不失一般性，这里以3次多项式为例，观察其d阶导数的多项式系数的变化规律：
  （1）3次多项式1阶导数的多项式系数相对于0阶导数的多项式系数，需要更新最后2个。最后一个系数乘以次数3，得到$3a_3$，倒数第二个系数乘以2(即次数减1)，得到$2a_2$。
  （2）3次多项式2阶导数的多项式系数相对于1阶导数的多项式系数，需要更新最后1个。最后一个系数乘以2(即次数减1)，得到$6a_3$。
  （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)