!************************************************************************* !* EVALUATE A POLYNOMIAL AND ITS DERIVATIVES BY HORNER'S METHOD * !* --------------------------------------------------------------------- * !* SAMPLE RUN: * !* * !* EVALUATE A POLYNOMIAL AND ITS DERIVATIVES BY HORNER'S METHOD * !* * !* Example: P(X) = X^4 + 2X^3 +3X^2 + 4X + 5 * !* * !* X P(X) P'(X) P"(X) P"'(X) * !* --------------------------------------------------------------------- * !* 0.00000E+00 0.50000E+01 0.40000E+01 0.60000E+01 0.12000E+02 * !* 0.20000E+00 0.59376E+01 0.54720E+01 0.88800E+01 0.16800E+02 * !* 0.40000E+00 0.72336E+01 0.76160E+01 0.12720E+02 0.21600E+02 * !* 0.60000E+00 0.90416E+01 0.10624E+02 0.17520E+02 0.26400E+02 * !* 0.80000E+00 0.11554E+02 0.14688E+02 0.23280E+02 0.31200E+02 * !* 0.10000E+01 0.15000E+02 0.20000E+02 0.30000E+02 0.36000E+02 * !* 0.12000E+01 0.19650E+02 0.26752E+02 0.37680E+02 0.40800E+02 * !* 0.14000E+01 0.25810E+02 0.35136E+02 0.46320E+02 0.45600E+02 * !* 0.16000E+01 0.33826E+02 0.45344E+02 0.55920E+02 0.50400E+02 * !* 0.18000E+01 0.44082E+02 0.57568E+02 0.66480E+02 0.55200E+02 * !* 0.20000E+01 0.57000E+02 0.72000E+02 0.78000E+02 0.60000E+02 * !* --------------------------------------------------------------------- * !* * !* Press any key to continue * !* * !* --------------------------------------------------------------------- * !* Reference: From Numath Library By Tuan Dang Trong in Fortran 77 * !* [BIBLI 18]. * !* * !* F90 Release 1.0 By J-P Moreau, Paris * !* (www.jpmoreau.fr) * !************************************************************************* PROGRAM TEST_HORNER REAL*8 A(25), RES(25) REAL*8 DX,X,X1,X2 LOGICAL FLG N = 4 !ORDER OF POLYNOMIAL NP = 11 !NUMBER OF POINTS !define coefficients A(1) = 1.d0 A(2) = 2.d0 A(3) = 3.d0 A(4) = 4.d0 A(5) = 5.d0 X1 = 0.d0 !X BEGIN X2 = 2.d0 !X END K=3 !NUMBER OF DERIVATIVES FLG=.TRUE. !EVALUATE ONLY DX = (X2-X1) / (NP-1) X = X1 - DX print *,' ' print *,' EVALUATE A POLYNOMIAL AND ITS DERIVATIVES BY HORNER''S METHOD' print *,' ' print *,' Example: P(X) = X^4 + 2X^3 +3X^2 + 4X + 5' print *,' ' print *,' X P(X) P''(X) P"(X) P"''(X)' print *,' ---------------------------------------------------------------------' DO I=1, NP X = X + DX CALL HORNER(N,A,X,K,FLG,RES) write(*,10) X, (RES(J), J=1,K+1) END DO print *,' ---------------------------------------------------------------------' print *,' ' 10 format(5E14.5) END SUBROUTINE HORNER (N,A,X0,K,B,R) !--------------------------------------------------------------------- ! EVALUATE A POLYNOMIAL AND ITS DERIVATIVES BY HORNER'S METHOD ! ! INPUTS: ! N ORDER OF POLYNOMIAL ! A VECTOR OF SIZE N+1 STORING THE COEFFICIENTS OF ! POLYNOMIAL IN DECREASING ORDERS OF POWERS ! I.E. P(X) = A(1)*X**N+A(2)*X**(N-1)+...+A(N)*X+A(N+1) ! X0 GIVEN ARGUMENT ! K MAXIMUM ORDER OF DERIVATIVES ASKED FOR ! B FLAG ! = .TRUE. EVALUATION ONLY ! = .FALSE. DETERMINATION OF POLYNOMIAL COEFFICIENTS ! NEAR X0 ! R VECTOR OF SIZE K+1 CONTAINS: ! - VALUES P(X0), P'(X0), P''(X0),..PK(X0), ! IF B IS TRUE AND IF K < N ! - THE N+1 COEFFICIENTS OF POLYNOMIAL ! Q(X-X0) = R(1)+R(2)*(X-X0)+...+R(N+1)*(X-X0)**(N), ! IF B IS FALSE AND IF K = N. ! ! REFERENCE: ! ALGORITHM 337, COLLECTED ALGORITHMS FROM CACM, W.PANKIEWICKZ !--------------------------------------------------------------------- REAL*8 A(*),R(*),X0 LOGICAL B RR = A(1) DO 10 I = 1,K+1 10 R(I) = RR DO 20 J = 2,N+1 R(1) = R(1)*X0+A(J) NMJ = N-J+1 IF (NMJ.GT.K) THEN L = K ELSE L = NMJ END IF DO 15 I = 2,L+1 R(I) = R(I)*X0+R(I-1) 15 CONTINUE 20 CONTINUE IF (B) THEN L = 1 DO 25 I = 2,K+1 L = L*(I-1) 25 R(I) = R(I)*L ELSE END IF RETURN END !end of file thorner.f90