Mercurial > hg > octave-nkf
changeset 5217:e88886a6934d
[project @ 2005-03-16 20:03:01 by jwe]
author | jwe |
---|---|
date | Wed, 16 Mar 2005 20:03:01 +0000 |
parents | 5ed60b8b1ac4 |
children | 4d036412ccca |
files | scripts/ChangeLog scripts/polynomial/polyder.m scripts/polynomial/polyderiv.m scripts/polynomial/polygcd.m |
diffstat | 4 files changed, 87 insertions(+), 77 deletions(-) [+] |
line wrap: on
line diff
--- a/scripts/ChangeLog +++ b/scripts/ChangeLog @@ -1,3 +1,9 @@ +2005-03-16 Paul Kienzle <pkienzle@users.sf.net> + + * polynomial/polyderiv.m : Accept a*b, a/b. Auto-reduce common terms. + * polynomial/polyder.m: Ditto. + * polynomial/polygcd.m: New function. + 2005-03-16 John W. Eaton <jwe@octave.org> * control/base/__stepimp__.m, control/base/bode.m,
--- a/scripts/polynomial/polyder.m +++ b/scripts/polynomial/polyder.m @@ -25,21 +25,19 @@ ## @end deftypefn ## Author: John W. Eaton -## Paul Kienzle <pkienzle@kienzle.powernet.co.uk> -## handle b/a and b*a function [q, r] = polyder (p, a) if (nargin == 1) q = polyderiv (p); - elseif (nargin==2) - if (nargout==2) - [q, r] = polyderiv (p,a); + elseif (nargin == 2) + if (nargout == 2) + [q, r] = polyderiv (p, a); else - q = polyderiv (p,a); + q = polyderiv (p, a); endif else - usage ("q=polyder(p) or q=polyder(b,a) or [q, r]=polyder(b,a)"); + usage ("q = polyder (p) or q = polyder (b, a) or [q, r] = polyder (b, a)"); endif endfunction
--- a/scripts/polynomial/polyderiv.m +++ b/scripts/polynomial/polyderiv.m @@ -33,63 +33,61 @@ ## Author: Tony Richardson <arichard@stark.cc.oh.us> ## Created: June 1994 ## Adapted-By: jwe -## Paul Kienzle <pkienzle@kienzle.powernet.co.uk> -## handle b/a and b*a function [q, r] = polyderiv (p, a) - if (nargin < 1 || nargin > 3) - usage ("q=polyderiv(p) or q=polyderiv(b,a) or [q, r]=polyderiv(b,a)"); - endif - - if (! isvector (p)) - error ("polyderiv: argument must be a vector"); - endif - - if (nargin == 2) - if (! isvector (a)) + if (nargin == 1 || nargin == 2) + if (! isvector (p)) error ("polyderiv: argument must be a vector"); endif - if (nargout == 1) - ## derivative of p*a returns a single polynomial - q = polyderiv(conv(p,a)); + if (nargin == 2) + if (! isvector (a)) + error ("polyderiv: argument must be a vector"); + endif + if (nargout == 1) + ## derivative of p*a returns a single polynomial + q = polyderiv (conv (p, a)); + else + ## derivative of p/a returns numerator and denominator + r = conv (a, a); + if (numel (p) == 1) + q = -p * polyderiv (a); + elseif (numel (a) == 1) + q = a * polyderiv (p); + else + q = conv (polyderiv (p), a) - conv (p, polyderiv (a)); + q = polyreduce (q); + endif + + ## remove common factors from numerator and denominator + x = polygcd (q, r); + if (length(x) != 1) + q = deconv (q, x); + r = deconv (r, x); + endif + + ## move all the gain into the numerator + q = q/r(1); + r = r/r(1); + endif else - ## derivative of p/a returns numerator and denominator - r = conv(a, a); - if numel(p) == 1 - q = -p * polyderiv(a); - elseif numel(a) == 1 - q = a * polyderiv(p); - else - q = conv(polyderiv(p),a) - conv(p,polyderiv(a)); - q = polyreduce(q); + lp = numel (p); + if (lp == 1) + q = 0; + return; + elseif (lp == 0) + q = []; + return; endif - ## remove common factors from numerator and denominator - x = polygcd(q,r); - if length(x)!=1 - q=deconv(q,x); - r=deconv(r,x); - endif + ## Force P to be a row vector. + p = p(:).'; - ## move all the gain into the numerator - q=q/r(1); - r=r/r(1); + q = p(1:(lp-1)) .* [(lp-1):-1:1]; endif else - lp = numel (p); - if (lp == 1) - q = 0; - return; - elseif (lp == 0) - q = []; - return; - end - - ## Force P to be a row vector. - p = p(:).'; - - q = p (1:(lp-1)) .* [(lp-1):-1:1]; + usage ("q = polyderiv (p) or q = polyderiv (b, a) or [q, r] = polyderiv (b, a)"); endif + endfunction
--- a/scripts/polynomial/polygcd.m +++ b/scripts/polynomial/polygcd.m @@ -37,29 +37,37 @@ ## @seealso{poly, polyinteg, polyderiv, polyreduce, roots, conv, deconv, ## residue, filter, polyval, and polyvalm} -function x = polygcd(b,a,tol) - if (nargin<2 || nargin>3) - usage("x=polygcd(b,a [,tol])"); - endif - if (nargin<3), tol=sqrt(eps); endif - if (length(a)==1 || length(b)==1) - if a==0, x=b; - elseif b==0, x=a; - else x=1; +function x = polygcd (b, a, tol) + + if (nargin == 2 || nargin == 3) + if (nargin == 2) + tol = sqrt (eps); endif - return; + if (length (a) == 1 || length (b) == 1) + if (a == 0) + x = b; + elseif (b == 0) + x = a; + else + x = 1; + endif + else + a /= a(1); + while (1) + [d, r] = deconv (b, a); + nz = find (abs (r) > tol); + if (isempty (nz)) + x = a; + break; + else + r = r(nz(1):length(r)); + endif + b = a; + a /= r(1); + endwhile + endif + else + usage ("x = polygcd (b, a [,tol])"); endif - a = a./a(1); - while (1) - [d, r] = deconv(b, a); - nz = find(abs(r)>tol); - if isempty(nz) - x = a; - return; - else - r = r(nz(1):length(r)); - endif - b = a; - a = r./r(1); - endwhile -endfunction \ No newline at end of file + +endfunction