1*DECK RPQR79
2      SUBROUTINE RPQR79 (NDEG, COEFF, ROOT, IERR, WORK)
3C***BEGIN PROLOGUE  RPQR79
4C***PURPOSE  Find the zeros of a polynomial with real coefficients.
5C***LIBRARY   SLATEC
6C***CATEGORY  F1A1A
7C***TYPE      SINGLE PRECISION (RPQR79-S, CPQR79-C)
8C***KEYWORDS  COMPLEX POLYNOMIAL, POLYNOMIAL ROOTS, POLYNOMIAL ZEROS
9C***AUTHOR  Vandevender, W. H., (SNLA)
10C***DESCRIPTION
11C
12C   Abstract
13C       This routine computes all zeros of a polynomial of degree NDEG
14C       with real coefficients by computing the eigenvalues of the
15C       companion matrix.
16C
17C   Description of Parameters
18C       The user must dimension all arrays appearing in the call list
19C            COEFF(NDEG+1), ROOT(NDEG), WORK(NDEG*(NDEG+2))
20C
21C    --Input--
22C      NDEG    degree of polynomial
23C
24C      COEFF   REAL coefficients in descending order.  i.e.,
25C              P(Z)= COEFF(1)*(Z**NDEG) + COEFF(NDEG)*Z + COEFF(NDEG+1)
26C
27C      WORK    REAL work array of dimension at least NDEG*(NDEG+2)
28C
29C   --Output--
30C      ROOT    COMPLEX vector of roots
31C
32C      IERR    Output Error Code
33C           - Normal Code
34C          0  means the roots were computed.
35C           - Abnormal Codes
36C          1  more than 30 QR iterations on some eigenvalue of the
37C             companion matrix
38C          2  COEFF(1)=0.0
39C          3  NDEG is invalid (less than or equal to 0)
40C
41C***REFERENCES  (NONE)
42C***ROUTINES CALLED  HQR, XERMSG
43C***REVISION HISTORY  (YYMMDD)
44C   800601  DATE WRITTEN
45C   890505  REVISION DATE from Version 3.2
46C   891214  Prologue converted to Version 4.0 format.  (BAB)
47C   900315  CALLs to XERROR changed to CALLs to XERMSG.  (THJ)
48C   911010  Code reworked and simplified.  (RWC and WRB)
49C***END PROLOGUE  RPQR79
50      REAL COEFF(*), WORK(*), SCALE
51      COMPLEX ROOT(*)
52      INTEGER NDEG, IERR, K, KH, KWR, KWI, KCOL
53C***FIRST EXECUTABLE STATEMENT  RPQR79
54      IERR = 0
55      IF (ABS(COEFF(1)) .EQ. 0.0) THEN
56         IERR = 2
57         CALL XERMSG ('SLATEC', 'RPQR79',
58     +      'LEADING COEFFICIENT IS ZERO.', 2, 1)
59         RETURN
60      ENDIF
61C
62      IF (NDEG .LE. 0) THEN
63         IERR = 3
64         CALL XERMSG ('SLATEC', 'RPQR79', 'DEGREE INVALID.', 3, 1)
65         RETURN
66      ENDIF
67C
68      IF (NDEG .EQ. 1) THEN
69         ROOT(1) = CMPLX(-COEFF(2)/COEFF(1),0.0)
70         RETURN
71      ENDIF
72C
73      SCALE = 1.0E0/COEFF(1)
74      KH = 1
75      KWR = KH+NDEG*NDEG
76      KWI = KWR+NDEG
77      KWEND = KWI+NDEG-1
78C
79      DO 10 K=1,KWEND
80         WORK(K) = 0.0E0
81   10 CONTINUE
82C
83      DO 20 K=1,NDEG
84         KCOL = (K-1)*NDEG+1
85         WORK(KCOL) = -COEFF(K+1)*SCALE
86         IF (K .NE. NDEG) WORK(KCOL+K) = 1.0E0
87   20 CONTINUE
88C
89      CALL HQR (NDEG,NDEG,1,NDEG,WORK(KH),WORK(KWR),WORK(KWI),IERR)
90C
91      IF (IERR .NE. 0) THEN
92         IERR = 1
93         CALL XERMSG ('SLATEC', 'CPQR79',
94     +      'NO CONVERGENCE IN 30 QR ITERATIONS.', 1, 1)
95         RETURN
96      ENDIF
97C
98      DO 30 K=1,NDEG
99         KM1 = K-1
100         ROOT(K) = CMPLX(WORK(KWR+KM1),WORK(KWI+KM1))
101   30 CONTINUE
102      RETURN
103      END
104