* E02RAF Example Program Text * Mark 16 Revised. NAG Copyright 1993. * .. Parameters .. INTEGER L, M, IA, IB, IC, IW PARAMETER (L=4,M=4,IA=L+1,IB=M+1,IC=IA+IB-1,IW=IB*(2*IB+3)) INTEGER NOUT PARAMETER (NOUT=6) LOGICAL SCALE PARAMETER (SCALE=.TRUE.) * .. Local Scalars .. INTEGER I, IFAIL * .. Local Arrays .. DOUBLE PRECISION AA(IA), BB(IB), CC(IC), DD(IA+IB), W(IW), + WORK(2*(L+M+1)), Z(2,L+M) * .. External Subroutines .. EXTERNAL C02AGF, E02RAF * .. Intrinsic Functions .. INTRINSIC DBLE * .. Executable Statements .. WRITE (NOUT,*) 'E02RAF Example Program Results' * Power series coefficients in CC CC(1) = 1.0D0 DO 20 I = 1, IC - 1 CC(I+1) = CC(I)/DBLE(I) 20 CONTINUE WRITE (NOUT,*) WRITE (NOUT,*) 'The given series coefficients are' WRITE (NOUT,99999) (CC(I),I=1,IC) IFAIL = 0 * CALL E02RAF(IA,IB,CC,IC,AA,BB,W,IW,IFAIL) * WRITE (NOUT,*) WRITE (NOUT,*) 'Numerator coefficients' WRITE (NOUT,99999) (AA(I),I=1,IA) WRITE (NOUT,*) WRITE (NOUT,*) 'Denominator coefficients' WRITE (NOUT,99999) (BB(I),I=1,IB) * Calculate zeros of the approximant using C02AGF * First need to reverse order of coefficients DO 40 I = 1, IA DD(IA-I+1) = AA(I) 40 CONTINUE IFAIL = 0 * CALL C02AGF(DD,L,SCALE,Z,WORK,IFAIL) * WRITE (NOUT,*) WRITE (NOUT,*) 'Zeros of approximant are at' WRITE (NOUT,*) WRITE (NOUT,*) ' Real part Imag part' WRITE (NOUT,99998) (Z(1,I),Z(2,I),I=1,L) * Calculate poles of the approximant using C02AGF * Reverse order of coefficients DO 60 I = 1, IB DD(IB-I+1) = BB(I) 60 CONTINUE IFAIL = 0 * CALL C02AGF(DD,M,SCALE,Z,WORK,IFAIL) * WRITE (NOUT,*) WRITE (NOUT,*) 'Poles of approximant are at' WRITE (NOUT,*) WRITE (NOUT,*) ' Real part Imag part' WRITE (NOUT,99998) (Z(1,I),Z(2,I),I=1,M) STOP * 99999 FORMAT (1X,5E13.4) 99998 FORMAT (1X,2E13.4) END