Code covered by the BSD License

# Polynomial division by convolution -- up to finite terms

### Feng Cheng Chang (view profile)

14 Apr 2008 (Updated )

Division of two polynomials by convolution to get up to K terms.

polydiv_Z(b,a,K,c)
```function [q,qc,c] = polydiv_Z(b,a,K,c)
%
%   Polynomial division by convolution.
%
%   Calculate inverse Z-transform (Polynomial division) up to K terms:
%          q(z) = b(z)/a(z),
%   where
%          b(z) = b(0) + ... + b(k)/z^k + ... + b(n)/z^n.
%          a(z) = a(0) + ... + a(k)/z^k + ... + a(m)/z^m.
%          q(z) = q(0) + ... + q(k)/z^k + ... + q(K)/z^K + ......
%
%   If coefficients of b(x) and a(x) are all integers, we set  c = 0,
%   so that the entire process involve mostly integer multiplications.
%   The round-off errors may thus be eliminated.
%
%   This code is similar to the code by Tamer Abdelazim Mellik,
%   "Calculate inverse Z-transform by long division."
%
%   By F C Chang  04/12/08  updated  07/07/11
%
if nargin < 4,  c = 1; end;
n = length(b);   m = length(a);
b = [b,zeros(1,K-1-n+m)];
if m == 1,  q = b/a;  qc = b;  c = a;    return;    end;
if c == 0,  w(1) = a(1)^K;  else  w(1) = 1; c = 1;  end;
for k = 2:K+1;
w(k) = [b(k-1),-a(min(k-1,m):-1:2)]*[w(1),w(max(2,k+1-m):k-1)]'/a(1);
end;
qc = w(2:K+1);
q = qc;
if c == 0, c = w(1), q = qc/c;  end;```