← Back to branch summary

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

« back to all changes in this revision

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