~ubuntu-branches/ubuntu/lucid/pdl/lucid

« back to all changes in this revision

Viewing changes to Lib/Slatec/slatec/dp1vlu.f

  • Committer: Bazaar Package Importer
  • Author(s): Ben Gertzfield
  • Date: 2002-04-08 18:47:16 UTC
  • Revision ID: james.westby@ubuntu.com-20020408184716-0hf64dc96kin3htp
Tags: upstream-2.3.2
ImportĀ upstreamĀ versionĀ 2.3.2

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
*DECK DP1VLU
 
2
      SUBROUTINE DP1VLU (L, NDER, X, YFIT, YP, A)
 
3
C***BEGIN PROLOGUE  DP1VLU
 
4
C***PURPOSE  Use the coefficients generated by DPOLFT to evaluate the
 
5
C            polynomial fit of degree L, along with the first NDER of
 
6
C            its derivatives, at a specified point.
 
7
C***LIBRARY   SLATEC
 
8
C***CATEGORY  K6
 
9
C***TYPE      DOUBLE PRECISION (PVALUE-S, DP1VLU-D)
 
10
C***KEYWORDS  CURVE FITTING, LEAST SQUARES, POLYNOMIAL APPROXIMATION
 
11
C***AUTHOR  Shampine, L. F., (SNLA)
 
12
C           Davenport, S. M., (SNLA)
 
13
C***DESCRIPTION
 
14
C
 
15
C     Abstract
 
16
C
 
17
C     The subroutine  DP1VLU  uses the coefficients generated by  DPOLFT
 
18
C     to evaluate the polynomial fit of degree  L , along with the first
 
19
C     NDER  of its derivatives, at a specified point.  Computationally
 
20
C     stable recurrence relations are used to perform this task.
 
21
C
 
22
C     The parameters for  DP1VLU  are
 
23
C
 
24
C     Input -- ALL TYPE REAL variables are DOUBLE PRECISION
 
25
C         L -      the degree of polynomial to be evaluated.  L  may be
 
26
C                  any non-negative integer which is less than or equal
 
27
C                  to  NDEG , the highest degree polynomial provided
 
28
C                  by  DPOLFT .
 
29
C         NDER -   the number of derivatives to be evaluated.  NDER
 
30
C                  may be 0 or any positive value.  If NDER is less
 
31
C                  than 0, it will be treated as 0.
 
32
C         X -      the argument at which the polynomial and its
 
33
C                  derivatives are to be evaluated.
 
34
C         A -      work and output array containing values from last
 
35
C                  call to  DPOLFT .
 
36
C
 
37
C     Output -- ALL TYPE REAL variables are DOUBLE PRECISION
 
38
C         YFIT -   value of the fitting polynomial of degree  L  at  X
 
39
C         YP -     array containing the first through  NDER  derivatives
 
40
C                  of the polynomial of degree  L .  YP  must be
 
41
C                  dimensioned at least  NDER  in the calling program.
 
42
C
 
43
C***REFERENCES  L. F. Shampine, S. M. Davenport and R. E. Huddleston,
 
44
C                 Curve fitting by polynomials in one variable, Report
 
45
C                 SLA-74-0270, Sandia Laboratories, June 1974.
 
46
C***ROUTINES CALLED  XERMSG
 
47
C***REVISION HISTORY  (YYMMDD)
 
48
C   740601  DATE WRITTEN
 
49
C   890531  Changed all specific intrinsics to generic.  (WRB)
 
50
C   890911  Removed unnecessary intrinsics.  (WRB)
 
51
C   891006  Cosmetic changes to prologue.  (WRB)
 
52
C   891006  REVISION DATE from Version 3.2
 
53
C   891214  Prologue converted to Version 4.0 format.  (BAB)
 
54
C   900315  CALLs to XERROR changed to CALLs to XERMSG.  (THJ)
 
55
C   900510  Convert XERRWV calls to XERMSG calls.  (RWC)
 
56
C   920501  Reformatted the REFERENCES section.  (WRB)
 
57
C***END PROLOGUE  DP1VLU
 
58
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
 
59
      INTEGER I,IC,ILO,IN,INP1,IUP,K1,K1I,K2,K3,K3P1,K3PN,K4,K4P1,K4PN,
 
60
     * KC,L,LM1,LP1,MAXORD,N,NDER,NDO,NDP1,NORD
 
61
      DOUBLE PRECISION A(*),CC,DIF,VAL,X,YFIT,YP(*)
 
62
      CHARACTER*8 XERN1, XERN2
 
63
C***FIRST EXECUTABLE STATEMENT  DP1VLU
 
64
      IF (L .LT. 0) GO TO 12
 
65
      NDO = MAX(NDER,0)
 
66
      NDO = MIN(NDO,L)
 
67
      MAXORD = A(1) + 0.5D0
 
68
      K1 = MAXORD + 1
 
69
      K2 = K1 + MAXORD
 
70
      K3 = K2 + MAXORD + 2
 
71
      NORD = A(K3) + 0.5D0
 
72
      IF (L .GT. NORD) GO TO 11
 
73
      K4 = K3 + L + 1
 
74
      IF (NDER .LT. 1) GO TO 2
 
75
      DO 1 I = 1,NDER
 
76
 1      YP(I) = 0.0D0
 
77
 2    IF (L .GE. 2) GO TO 4
 
78
      IF (L .EQ. 1) GO TO 3
 
79
C
 
80
C L IS 0
 
81
C
 
82
      VAL = A(K2+1)
 
83
      GO TO 10
 
84
C
 
85
C L IS 1
 
86
C
 
87
 3    CC = A(K2+2)
 
88
      VAL = A(K2+1) + (X-A(2))*CC
 
89
      IF (NDER .GE. 1) YP(1) = CC
 
90
      GO TO 10
 
91
C
 
92
C L IS GREATER THAN 1
 
93
C
 
94
 4    NDP1 = NDO + 1
 
95
      K3P1 = K3 + 1
 
96
      K4P1 = K4 + 1
 
97
      LP1 = L + 1
 
98
      LM1 = L - 1
 
99
      ILO = K3 + 3
 
100
      IUP = K4 + NDP1
 
101
      DO 5 I = ILO,IUP
 
102
 5      A(I) = 0.0D0
 
103
      DIF = X - A(LP1)
 
104
      KC = K2 + LP1
 
105
      A(K4P1) = A(KC)
 
106
      A(K3P1) = A(KC-1) + DIF*A(K4P1)
 
107
      A(K3+2) = A(K4P1)
 
108
C
 
109
C EVALUATE RECURRENCE RELATIONS FOR FUNCTION VALUE AND DERIVATIVES
 
110
C
 
111
      DO 9 I = 1,LM1
 
112
        IN = L - I
 
113
        INP1 = IN + 1
 
114
        K1I = K1 + INP1
 
115
        IC = K2 + IN
 
116
        DIF = X - A(INP1)
 
117
        VAL = A(IC) + DIF*A(K3P1) - A(K1I)*A(K4P1)
 
118
        IF (NDO .LE. 0) GO TO 8
 
119
        DO 6 N = 1,NDO
 
120
          K3PN = K3P1 + N
 
121
          K4PN = K4P1 + N
 
122
 6        YP(N) = DIF*A(K3PN) + N*A(K3PN-1) - A(K1I)*A(K4PN)
 
123
C
 
124
C SAVE VALUES NEEDED FOR NEXT EVALUATION OF RECURRENCE RELATIONS
 
125
C
 
126
        DO 7 N = 1,NDO
 
127
          K3PN = K3P1 + N
 
128
          K4PN = K4P1 + N
 
129
          A(K4PN) = A(K3PN)
 
130
 7        A(K3PN) = YP(N)
 
131
 8      A(K4P1) = A(K3P1)
 
132
 9      A(K3P1) = VAL
 
133
C
 
134
C NORMAL RETURN OR ABORT DUE TO ERROR
 
135
C
 
136
 10   YFIT = VAL
 
137
      RETURN
 
138
C
 
139
   11 WRITE (XERN1, '(I8)') L
 
140
      WRITE (XERN2, '(I8)') NORD
 
141
      CALL XERMSG ('SLATEC', 'DP1VLU',
 
142
     *   'THE ORDER OF POLYNOMIAL EVALUATION, L = ' // XERN1 //
 
143
     *   ' REQUESTED EXCEEDS THE HIGHEST ORDER FIT, NORD = ' // XERN2 //
 
144
     *   ', COMPUTED BY DPOLFT -- EXECUTION TERMINATED.', 8, 2)
 
145
      RETURN
 
146
C
 
147
   12 CALL XERMSG ('SLATEC', 'DP1VLU',
 
148
     +   'INVALID INPUT PARAMETER.  ORDER OF POLYNOMIAL EVALUATION ' //
 
149
     +   'REQUESTED IS NEGATIVE.', 2, 2)
 
150
      RETURN
 
151
      END