|
559
|
1 function r = roots(c) |
|
|
2 #Find the roots of a polynomial. |
|
|
3 # |
|
|
4 #In octave, a polynomial is represented by it's coefficients (arranged |
|
|
5 #in descending order). For example, a vector c of length n+1 corresponds |
|
|
6 #to the following nth order polynomial |
|
|
7 # |
|
|
8 # p(x) = c(1) x^n + ... + c(n) x + c(n+1). |
|
|
9 # |
|
|
10 #roots(c) will return a vector r of length n such that |
|
|
11 # |
|
|
12 # p(x) = c(1) [ (x - r(1)) * (x - r(2)) * ... * (x - r(n)) ] |
|
|
13 # |
|
|
14 #roots and poly are inverse functions to within a scaling factor. |
|
|
15 # |
|
|
16 #SEE ALSO: poly, roots, conv, deconv, residue, filter, polyval, polyvalm, |
|
|
17 # polyderiv, polyinteg |
|
|
18 |
|
|
19 # Author: |
|
|
20 # Tony Richardson |
|
|
21 # amr@mpl.ucsd.edu |
|
|
22 # June 1994 |
|
|
23 |
|
|
24 if(nargin != 1) |
|
|
25 error("usage: roots(vector)"); |
|
|
26 endif |
|
|
27 |
|
|
28 if(is_matrix(c)) |
|
|
29 error("argument must be a vector."); |
|
|
30 endif |
|
|
31 |
|
|
32 n = length(c); |
|
|
33 |
|
|
34 if(is_scalar(c) || n == 0) |
|
|
35 r = []; |
|
|
36 return; |
|
|
37 endif |
|
|
38 |
|
|
39 # Ensure that c is a row vector. |
|
|
40 if(rows(c) > 1) |
|
|
41 c = reshape(c,1,n); |
|
|
42 endif |
|
|
43 |
|
|
44 # We could replace this with a call to compan, but it's faster to |
|
|
45 # just reproduce the code here. |
|
|
46 A = diag(ones(n-2,1),-1); |
|
|
47 A(1,:) = -c(2:n)/c(1); |
|
|
48 |
|
|
49 r = eig(A); |
|
|
50 |
|
|
51 # Sort roots in order by decreasing magnitude. |
|
|
52 [mr i] = sort(abs(r)); |
|
|
53 r = r(i(length(i):-1:1)); |
|
|
54 |
|
|
55 endfunction |
