pchdf.f 3.1 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106
  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