dchfie.f 3.3 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109
  1. *DECK DCHFIE
  2. DOUBLE PRECISION FUNCTION DCHFIE (X1, X2, F1, F2, D1, D2, A, B)
  3. C***BEGIN PROLOGUE DCHFIE
  4. C***SUBSIDIARY
  5. C***PURPOSE Evaluates integral of a single cubic for DPCHIA
  6. C***LIBRARY SLATEC (PCHIP)
  7. C***TYPE DOUBLE PRECISION (CHFIE-S, DCHFIE-D)
  8. C***AUTHOR Fritsch, F. N., (LLNL)
  9. C***DESCRIPTION
  10. C
  11. C DCHFIE: Cubic Hermite Function Integral Evaluator.
  12. C
  13. C Called by DPCHIA to evaluate the integral of a single cubic (in
  14. C Hermite form) over an arbitrary interval (A,B).
  15. C
  16. C ----------------------------------------------------------------------
  17. C
  18. C Calling sequence:
  19. C
  20. C DOUBLE PRECISION X1, X2, F1, F2, D1, D2, A, B
  21. C DOUBLE PRECISION VALUE, DCHFIE
  22. C
  23. C VALUE = DCHFIE (X1, X2, F1, F2, D1, D2, A, B)
  24. C
  25. C Parameters:
  26. C
  27. C VALUE -- (output) value of the requested integral.
  28. C
  29. C X1,X2 -- (input) endpoints if interval of definition of cubic.
  30. C
  31. C F1,F2 -- (input) function values at the ends of the interval.
  32. C
  33. C D1,D2 -- (input) derivative values at the ends of the interval.
  34. C
  35. C A,B -- (input) endpoints of interval of integration.
  36. C
  37. C***SEE ALSO DPCHIA
  38. C***ROUTINES CALLED (NONE)
  39. C***REVISION HISTORY (YYMMDD)
  40. C 820730 DATE WRITTEN
  41. C 820805 Converted to SLATEC library version.
  42. C 870707 Corrected subroutine name from DCHIV to DCHFIV.
  43. C 870813 Minor cosmetic changes.
  44. C 890411 1. Added SAVE statements (Vers. 3.2).
  45. C 2. Added SIX to DOUBLE PRECISION declaration.
  46. C 890411 REVISION DATE from Version 3.2
  47. C 891214 Prologue converted to Version 4.0 format. (BAB)
  48. C 900315 CALLs to XERROR changed to CALLs to XERMSG. (THJ)
  49. C 900328 Added TYPE section. (WRB)
  50. C 910408 Updated AUTHOR section in prologue. (WRB)
  51. C 930503 Corrected to set VALUE=0 when IERR.ne.0. (FNF)
  52. C 930504 Eliminated IERR and changed name DCHFIV to DCHFIE. (FNF)
  53. C***END PROLOGUE DCHFIE
  54. C
  55. C Programming notes:
  56. C 1. There is no error return from this routine because zero is
  57. C indeed the mathematically correct answer when X1.EQ.X2 .
  58. C**End
  59. C
  60. C DECLARE ARGUMENTS.
  61. C
  62. DOUBLE PRECISION X1, X2, F1, F2, D1, D2, A, B
  63. C
  64. C DECLARE LOCAL VARIABLES.
  65. C
  66. DOUBLE PRECISION DTERM, FOUR, FTERM, H, HALF, PHIA1, PHIA2,
  67. * PHIB1, PHIB2, PSIA1, PSIA2, PSIB1, PSIB2, SIX, TA1, TA2,
  68. * TB1, TB2, THREE, TWO, UA1, UA2, UB1, UB2
  69. SAVE HALF, TWO, THREE, FOUR, SIX
  70. C
  71. C INITIALIZE.
  72. C
  73. DATA HALF/.5D0/, TWO/2.D0/, THREE/3.D0/, FOUR/4.D0/, SIX/6.D0/
  74. C
  75. C VALIDITY CHECK INPUT.
  76. C
  77. C***FIRST EXECUTABLE STATEMENT DCHFIE
  78. IF (X1 .EQ. X2) THEN
  79. DCHFIE = 0
  80. ELSE
  81. H = X2 - X1
  82. TA1 = (A - X1) / H
  83. TA2 = (X2 - A) / H
  84. TB1 = (B - X1) / H
  85. TB2 = (X2 - B) / H
  86. C
  87. UA1 = TA1**3
  88. PHIA1 = UA1 * (TWO - TA1)
  89. PSIA1 = UA1 * (THREE*TA1 - FOUR)
  90. UA2 = TA2**3
  91. PHIA2 = UA2 * (TWO - TA2)
  92. PSIA2 = -UA2 * (THREE*TA2 - FOUR)
  93. C
  94. UB1 = TB1**3
  95. PHIB1 = UB1 * (TWO - TB1)
  96. PSIB1 = UB1 * (THREE*TB1 - FOUR)
  97. UB2 = TB2**3
  98. PHIB2 = UB2 * (TWO - TB2)
  99. PSIB2 = -UB2 * (THREE*TB2 - FOUR)
  100. C
  101. FTERM = F1*(PHIA2 - PHIB2) + F2*(PHIB1 - PHIA1)
  102. DTERM = ( D1*(PSIA2 - PSIB2) + D2*(PSIB1 - PSIA1) )*(H/SIX)
  103. C
  104. DCHFIE = (HALF*H) * (FTERM + DTERM)
  105. ENDIF
  106. C
  107. RETURN
  108. C------------- LAST LINE OF DCHFIE FOLLOWS -----------------------------
  109. END