1 | *DECK PCHDF |
---|
2 | REAL FUNCTION PCHDF (K, X, S, IERR) |
---|
3 | C***BEGIN PROLOGUE PCHDF |
---|
4 | C***SUBSIDIARY |
---|
5 | C***PURPOSE Computes divided differences for PCHCE and PCHSP |
---|
6 | C***LIBRARY SLATEC (PCHIP) |
---|
7 | C***TYPE SINGLE PRECISION (PCHDF-S, DPCHDF-D) |
---|
8 | C***AUTHOR Fritsch, F. N., (LLNL) |
---|
9 | C***DESCRIPTION |
---|
10 | C |
---|
11 | C PCHDF: PCHIP Finite Difference Formula |
---|
12 | C |
---|
13 | C Uses a divided difference formulation to compute a K-point approx- |
---|
14 | C imation to the derivative at X(K) based on the data in X and S. |
---|
15 | C |
---|
16 | C Called by PCHCE and PCHSP to compute 3- and 4-point boundary |
---|
17 | C derivative approximations. |
---|
18 | C |
---|
19 | C ---------------------------------------------------------------------- |
---|
20 | C |
---|
21 | C On input: |
---|
22 | C K is the order of the desired derivative approximation. |
---|
23 | C K must be at least 3 (error return if not). |
---|
24 | C X contains the K values of the independent variable. |
---|
25 | C X need not be ordered, but the values **MUST** be |
---|
26 | C distinct. (Not checked here.) |
---|
27 | C S contains the associated slope values: |
---|
28 | C S(I) = (F(I+1)-F(I))/(X(I+1)-X(I)), I=1(1)K-1. |
---|
29 | C (Note that S need only be of length K-1.) |
---|
30 | C |
---|
31 | C On return: |
---|
32 | C S will be destroyed. |
---|
33 | C IERR will be set to -1 if K.LT.2 . |
---|
34 | C PCHDF will be set to the desired derivative approximation if |
---|
35 | C IERR=0 or to zero if IERR=-1. |
---|
36 | C |
---|
37 | C ---------------------------------------------------------------------- |
---|
38 | C |
---|
39 | C***SEE ALSO PCHCE, PCHSP |
---|
40 | C***REFERENCES Carl de Boor, A Practical Guide to Splines, Springer- |
---|
41 | C Verlag, New York, 1978, pp. 10-16. |
---|
42 | C***ROUTINES CALLED XERMSG |
---|
43 | C***REVISION HISTORY (YYMMDD) |
---|
44 | C 820503 DATE WRITTEN |
---|
45 | C 820805 Converted to SLATEC library version. |
---|
46 | C 870813 Minor cosmetic changes. |
---|
47 | C 890411 Added SAVE statements (Vers. 3.2). |
---|
48 | C 890411 REVISION DATE from Version 3.2 |
---|
49 | C 891214 Prologue converted to Version 4.0 format. (BAB) |
---|
50 | C 900315 CALLs to XERROR changed to CALLs to XERMSG. (THJ) |
---|
51 | C 900328 Added TYPE section. (WRB) |
---|
52 | C 910408 Updated AUTHOR and DATE WRITTEN sections in prologue. (WRB) |
---|
53 | C 920429 Revised format and order of references. (WRB,FNF) |
---|
54 | C 930503 Improved purpose. (FNF) |
---|
55 | C***END PROLOGUE PCHDF |
---|
56 | C |
---|
57 | C**End |
---|
58 | C |
---|
59 | C DECLARE ARGUMENTS. |
---|
60 | C |
---|
61 | INTEGER K, IERR |
---|
62 | REAL X(K), S(K) |
---|
63 | C |
---|
64 | C DECLARE LOCAL VARIABLES. |
---|
65 | C |
---|
66 | INTEGER I, J |
---|
67 | REAL VALUE, ZERO |
---|
68 | SAVE ZERO |
---|
69 | DATA ZERO /0./ |
---|
70 | C |
---|
71 | C CHECK FOR LEGAL VALUE OF K. |
---|
72 | C |
---|
73 | C***FIRST EXECUTABLE STATEMENT PCHDF |
---|
74 | IF (K .LT. 3) GO TO 5001 |
---|
75 | C |
---|
76 | C COMPUTE COEFFICIENTS OF INTERPOLATING POLYNOMIAL. |
---|
77 | C |
---|
78 | DO 10 J = 2, K-1 |
---|
79 | DO 9 I = 1, K-J |
---|
80 | S(I) = (S(I+1)-S(I))/(X(I+J)-X(I)) |
---|
81 | 9 CONTINUE |
---|
82 | 10 CONTINUE |
---|
83 | C |
---|
84 | C EVALUATE DERIVATIVE AT X(K). |
---|
85 | C |
---|
86 | VALUE = S(1) |
---|
87 | DO 20 I = 2, K-1 |
---|
88 | VALUE = S(I) + VALUE*(X(K)-X(I)) |
---|
89 | 20 CONTINUE |
---|
90 | C |
---|
91 | C NORMAL RETURN. |
---|
92 | C |
---|
93 | IERR = 0 |
---|
94 | PCHDF = VALUE |
---|
95 | RETURN |
---|
96 | C |
---|
97 | C ERROR RETURN. |
---|
98 | C |
---|
99 | 5001 CONTINUE |
---|
100 | C K.LT.3 RETURN. |
---|
101 | IERR = -1 |
---|
102 | CALL XERMSG ('SLATEC', 'PCHDF', 'K LESS THAN THREE', IERR, 1) |
---|
103 | PCHDF = ZERO |
---|
104 | RETURN |
---|
105 | C------------- LAST LINE OF PCHDF FOLLOWS ------------------------------ |
---|
106 | END |
---|